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Dynamics  of  a  system  of  many  bodies  in  space  is  formulated  in  a  Hamil¬ 
tonian  setting.  Typically  there  are  symmetry  groups  associated  to  such  problems, 
and  one  can  reduce  the  phase  space  by  these  symmetry  groups.  Dynamics  in  the 
reduced  phase  space  is  determined  by  appropriate  Poisson  structure.  Equilibria 
for  specific  cases  are  obtained  and  their  stability  examined  using  energy-Casimir 
method.  Nonlinear  control  techniques  including  exact  linearization  are  applied 
successfully.  A  global  controllability  theorem  is  proved  and  feedback  stabilization 
is  studied.  Applications  to  robotics  and  multibody  systems  in  space  are  discussed. 
Symbolic  computational  tools  were  used  extensively  in  the  research.  OOPSS  - 
an  Object  Oriented  Planar  System  Simulator  package  has  been  developed  to  au¬ 
tomatically  generate,  analyze,  simulate  and  graphically  display  the  dynamics  of 
planar  multibody  systems. 
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CHAPTER 

_ ONE 

INTRODUCTION 


With  the  advent  of  large  sophisticated  spacecraft  with  multiple  missions,  the  num¬ 
ber  and  size  of  the  attachments  on  the  main  satellite  body  has  grown.  Control  of  such 
multibody  spacecraft  is  difficult  since  the  kinematics  and  the  dynamics  in  such  systems 
are  coupled.  A  better  and  deeper  understanding  of  the  dynamics  is  imperative  for 
pragmatic  design  of  such  multi-body  systems  in  space.  Equipped  with  this  knowledge, 
we  would  be  able  to  design  control  systems  to  realize  the  mission  objectives. 

Spacecraft  in  general  are  collections  of  rigid  and  elastic  bodies,  sensors  and  actua¬ 
tors.  Complicating  issues  such  as  unmodeled/unmodelable  dynamics,  nonlinearities, 
reconfiguration,  external  disturbances,  material  properties  (including  material  aging 
effects)  are  impossible  to  treat  in  a  completely  rigorous  and  deterministic  way  [29].  In¬ 
herent  inability  to  model  a  system  perfectly  has  motivated  feedback  control.  However, 
a  significant  degree  of  modeling  is  absolutely  necessary  to  develop  such  control. 

Currently  large  antennae,  instrumentation,  solar  arrays  and  various  manipulator 
systems  are  being  designed  as  possible  attachments.  Examples  include  SASP  (Science 
and  Application  Systems  Package),  LDF  (Large  Deflector  Facility),  TWS  (Tele-robot 
Work  System)  etc..  These  attachments  are  usually  connected  to  the  main  body  by 
revolute  joints.  The  movement  of  the  one  body  relative  to  another  induces  reaction 
torques  and  forces  on  the  main  satellite  body  and  also  on  the  connected  bodies. 
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Extensive  use  of  remote  manipulators  and  robots  for  assembly,  construction,  re¬ 
pair,  and  servicing  of  satellites  is  being  considered  seriously  [28],  [6], [61].  Recent 
developments  in  automation  and  robotics  have  increased  the  importance  of  their  ap¬ 
plications  to  future  space  endeavors.  It  has  already  been  demonstrated  that  the  RMS 
(  Remote  Manipulator  System)  on  the  space  shuttle  could  be  used  to  capture  and 
service  a  satellite  [7].  But  suitable  corrections  have  to  be  made  on  the  Space  Shuttle 
to  compensate  for  the  arm  movements.  The  MRMS  (Mobile  RMS)  on  the  proposed 
Space  Station  will  be  used  for  assembly  and  construction  of  Space  Station  itself  [62]. 
The  reaction  torques  and  forces  on  the  Space  Station  due  to  the  MRMS  movement 
may  cause  the  center  of  mass  of  the  station  to  translate  and  may  rotate  the  station 
about  its  center  of  mass  [39].  Tele-operation  of  manipulator  arms  is  expected  to  be 
crucial  in  replacing  a  sizable  amount  of  EVA  (Extra  Vehicular  Activity).  The  ROBIN 
[6]  and  the  Grumman  tele-robot  [61]  are  typical  Telerobotic  Work  System  (TWS) 
configurations.  These  missions  involve  complex  tasks  and  call  for  a  higher  degree  of 
manipulation  and  control  than  is  available  now.  Furthermore,  the  mission  objectives 
lay  down  stringent  specifications  on  the  attitude  and  orbit  control  of  such  spacecraft. 

The  above  sets  of  complex  tasks  and  and  requirements  have  generated  a  number  of 
challenging  problems  in  the  area  of  dynamics  and  control  of  multi-body  spacecraft.  A 
thorough  understanding  of  the  dynamics  of  such  systems  is  required  for  good  control 
system  design.  Knowledge  of  the  relative  equilibria  and  their  stability  is  essential. 
Movement  of  these  equilibria  under  feed-back  control  needs  to  be  examined  closely. 
Thus  a  need  for  a  complete  understanding  of  such  dynamic  systems  is  identified. 

With  the  growing  importance  of  spacecraft  configurations  with  many  bodies  con¬ 
nected  together  by  revolute  joints,  problems  involving  re-orientation  of  such  space¬ 
craft  wherein  the  bodies  are  comparable  in  terms  of  their  physical  parameters  have 
become  important  [46].  The  spacecraft  slewing  (large  angle)  maneuver  problem  is 
inherently  nonlinear  owing  to  nonlinear  kinematics  and  gyroscopic  coupling  effects 
[29].  Time-optimal  slewing  of  spacecraft  has  been  considered  in  the  Lagrangian  and 
Newton-Euler  setting  and  the  solutions  provided  range  from  using  numerical  meth¬ 
ods  [55] (like  relaxation  and  continuation  methods),  to  using  feedback  linearization 
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[23]  and  then  solving  for  the  control  [11]. 

It  is  the  purpose  of  this  thesis  to  gain  better  understanding  and  insight  into  the 
dynamics  and  control  of  multi-body  systems.  We  deal  with  planar  multibody  systems 
connected  in  the  form  of  open  kinematic  chain  (no  closed  loops).  Individual  bodies 
are  modeled  as  rigid  bodies.  The  thesis  is  organized  as  follows  :  Chapter  two  gives 
the  basic  geometric  results  which  facilitate  understanding  of  the  modeling  aspects. 
Chapter  three  deals  with  the  construction  of  dynamical  models  of  planar  multibody 
systems  in  Hamiltonian  setting.  In  Chapter  four  we  give  some  examples  of  planar 
multibody  systems,  -  two-body,  three-body  and  chain-body  systems  in  space  ,  and 
planar  body  systems  on  ground.  We  examine  the  equilibria  of  these  dynamical  systems 
and  their  stability  using  the  energy- Casimir  method  in  Chapter  five.  Chapter  six 
deals  with  controllability,  feedback  linearization  and  feedback  stabilization  issues.  In 
Chapter  seven  we  present  some  symbolic  and  graphic  tools  we  developed  as  a  part  of 
our  research,  in  particular  OOPSS  (Object  Oriented  Planar  System  Simulator). 

Theorems,  Lemmas,  Corollaries,  Definitions  and  Remarks  are  numbered  individ¬ 
ually,  with  the  first  number  representing  the  chapter  they  occur  and  the  second  the 
section.  For  example,  Corollary  3.5.1  is  a  Corollary  in  Chapter  3,  Section  5  and  is 
the  first  Corollary  in  that  section.  The  next  Corollary  in  the  same  section  will  be 
Corollary  3.5.2. 
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CHAPTER 

_ TWO 

GEOMETRIC  PRELIMINARIES 


In  this  chapter  we  discuss  some  mathematical  definitions  and  concepts  as  related 
to  the  subject  material  on  hand.  We  assume  that  the  reader  is  familiar  with  the 
basic  concepts  such  as  manifolds,  vector  fields  (see  Abraham,  Marsden  and  Ratiu  [3] 
Chapter  3  and  4  for  a  good  exposition  of  these  basic  concepts)  and  symmetries  (see 
Abraham  and  Marsden  [2],  Chapter  4,  p.  341).  As  we  proceed  in  this  chapter  we  give 
necessary  ideas  and  definitions  or  refer  the  reader  to  other  sources. 

Let  M  be  a  smooth  manifold.  Let  C°°(U )  be  the  set  of  all  real  valued  smooth 
functions  in  a  neighborhood  U  of  p,  where  U  C  M.  There  exists  a  tangent  vector 
at  p  to  every  smooth  curve  in  the  neighborhood  U  of  p  and  passing  through  p.  The 
tangent  space  TPM  to  M  at  p  is  the  set  of  all  tangent  vectors  to  all  the  smooth 
curves  in  the  neighborhood  U  of  p  and  passing  through  p.  The  tangent  bundle  TM 
is  defined  on  M  to  be 


TM  =  {TPM  :  Vp  G  M}  , 

which  is  the  union  of  tangent  spaces  of  all  points  on  M.  The  linear  dual  of  tangent 
space  TPM  is  the  cotangent  space  T*M  composed  of  covectors  at  p.  Similarly  the 
dual  of  the  tangent  bundle  TM  is  the  cotangent  bundle  T*M. 
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2.1  Lie  Groups 

The  following  material  has  been  developed  as  in  Abraham  and  Marsden  [2]  Chap¬ 
ter  4. 

A  Lit  group  G  is  a  smooth  manifold  with  a  group  structure  and  the  group  oper¬ 
ations  of  multiplication  and  inversion  being  smooth.  The  vector  space  TtG  with  the 
Lit  algtbra  structure 

le,»7]  =  [*«,*„]  for  t,rjeTeG 

is  called  the  Lit  algtbra  of  G  and  is  denoted  by  g.  Here  X $  is  the  vector  field  obtained 
by  translating  £  to  all  of  G,  and  t  is  the  identity  of  G. 

Action  of  a  Group  :  An  action  of  a  Lit  Group  G  on  a  smooth  manifold  Af  is  a 
smooth  mapping  such  that  Vx  £  Af, 

(i)  $(e,  x)  =  x, 

(ii)  for  every  g,h  £  G,  4>[g,<t>{h,x))  —  if>(gh,x). 

The  action  <bg  of  g  £  G,  is  defined  as  :  Af  — ►  Af  :  x  *-*■  $(g,  x).  The  orbit  or 

flow  of  x  is  given  by  G  •  x  =  {$,(x)  :  g  £  G}. 

Infinitesimal  Generator  :  If  £  £  TeG ,  then  :  R  x  Af  — ►  Af  :  (t,x)  (->• 

^(exp^^jx)1  is  an  II  -  action  on  M  i.e.,  is  a  flow  on  Af.  The  corresponding 

■> 

vector  field  on  Af  given  by 

&/(*)  =  -^$(exp(t£),z) 

ai  t= o 

is  called  the  infinitesimal  generator  of  the  action  corresponding  to  £. 

Pull  back  &  Push  forward  :  Given  two  manifolds  Af  and  N  with  7[N)  = 
{set  of  all  real  valued  C°°  functions  on N}  and  X  (Af)  =  {vector  fields  on  M}  and  <f>  : 
Af  — ►  N  is  a  diffeomorphism,  the  pull  back  of  /  by  tf>  is 

4>'f  = 

1see  Abraham  and  Marsden  [2]  page  256  for  the  definition  of  exponential  mapping. 
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The  push  forward  of  X  by  <f>  is  defined  as 


<f>,x  =  T4>  oXof1  e  X(N). 


2.2  Mechanics  on  Manifolds 


The  following  material  is  patterned  after  Omohundro  [43]  and  Abraham  and  Mars- 
den  [2]  Chapter  3  and  4. 

2.2.1  Dynamics  on  Manifolds 

The  evolution  of  mechanical  systems  could  be  described  classically  in  terms  of 
the  generalized  coordinates  ft  coordinatizing  the  configuration  space  Q  ,  and  a  scalar 
quantity  L  called  the  Lagrangian  defined  on  the  tangent  space  TQ  of  Q  i.e., 


L  :  TQ-+1R.  :  (ft,  ft)  ^  £(ft,ft). 


The  Lagrangian  is  the  difference  between  the  kinetic  and  the  potential  energies  of 
the  system.  The  dynamics  of  the  system  is  given  by 


=  0. 


The  evolution  of  classical  mechanical  systems  could  also  be  described  in  terms 
of  the  generalized  coordinates  ft,  the  conjugate  momenta  p,-  and  the  corresponding 
Hamiltonian  H  which  maps  the  cotangent  space  into  the  real  line  i.e., 


H  :  T*Q  -+  R  :  (ft,p<)  H{qi,pi) 


The  dynamics  of  the  system  could  be  represented  as 


ft 


dH 
dPi' 
dH 
"  dqx  ' 


V 
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In  association  with  the  Hamiltonian  description,  one  has  the  notion  of  a  Poisson 
bracket  {•,  •}.  The  Poisson  bracket 


at 


{ft  9}  =  X) 
»=1 


df  dg  df  dg 


dqidpi  dpi  dq{  \ 

of  two  functions  of  qi  and  p,-,  defines  the  Poisson  structure  on  the  manifold  T'Q.  The 
above  bracket  is  known  as  the  canonical  bracket. 

A  Poisson  bracket  is  a  bilinear  map  from  pairs  of  functions  to  functions  which 
makes  the  space  of  smooth  functions  on  T*Q  into  a  Lie  Algebra,  and  acts  on 
products  like  a  derivative  does.  The  axioms  obeyed  by  the  Poisson  Brackets  are  : 


Bilinearity:  {afx  +  bf2,g}  =  a  {fx,g}  +  b{f2,g}  , 

Anti-symmetry  :  {ft  9}  =  —  {9tf}  » 

Jacobi’s  identity  :  {ft  {9th}}  +  {g,  {h,  f}}  +  {h,{f,g}}  =  0, 

Derivation  property  :  {f,gh}  =  {f,g}h  +  {f,h}g. 

A  manifold  P  with  a  Poisson  structure  is  a  Poisson  manifold ',  on  which  the 
Hamiltonian  is  a  function.  The  dynamics  evolve  on  this  Poisson  manifold.  The 
evolution  (dynamical)  equations  could  also  be  represented  in  terms  of  the  bracket  as 


q{  =  { qi,H }, 

Pi  =  {Pi,H}. 

Any  function  which  commutes  with  every  function  on  the  Poisson  manifold  is  auto¬ 
matically  a  constant  of  motion  and  is  called  a  Casimir  function.  Thus  <j)  is  a  Casimir 
function  iff, 


{<f>,g}  =  0 

If  there  exists  a  closed  nondegenerate  two-form  *,  w  called  a  symplectic  structure 
then  the  Poisson  manifold  P  is  known  as  the  symplectic  manifold  (P,  w).  Here 
u  :  Vect (P)  x  Vect(P)  — *■  P(P)  is  a  bilinear  skew-symmetric  map ,  F(P)  is  the 
space  of  smooth  real  valued  function  on  P  and  Vect(P)  is  the  space  of  smooth  vector 
fields  on  P. 

2If  w[X,  Y)  =  0  V  Y  6  Vect(P)  implies  X  =  0,  then  w  is  nondegenerate. 
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Let  6  F(P)  an<3  let  X /  be  &  vector  field  of  F  such  that 

u{Xf,Y)  =  df(Y)  VY  €  Vect(P) 

Now  the  Poisson  bracket  could  be  defined  as 

F(P)  x  F(P)  -  F(P) :  (/,})  >->  {/,s}  =  u[Xr,X,) 

From  the  properties  of  w  it  follows  that  (F(P),  {•,  •})  is  a  Lie  algebra.  Since  u  is  a 
nondegenerate  form,  the  associated  Poisson  structure  is  said  to  be  nonsingular.  So 
the  concept  of  Poisson  manifold  generalizes  the  notion  of  symplectic  manifold. 

A  natural  symplectic  structure  arises  from  Lagrangian  mechanical  systems  on  the 
configuration  space  Q.  The  Lagrangian  description  is  based  on  the  tangent  bundle 
TQ.  The  Hamiltonian  lives  on  the  cotangent  bundle  T*Q  (momenta  being  deriva¬ 
tives  of  L  with  respect  to  velocity  are  naturally  dual  to  velocities  and  thus  are  cov¬ 
ectors).  The  fiber  derivative  FL  relates  points  on  the  tangent  bundle  to  points 
on  the  cotangent  bundle.  The  following  diagram  helps  keep  the  relations  in  mind. 

FL  L 

T*Q  ◄ -  TQ - 


Here  E  is  called  the  energy  of  the  Lagrangian.  T*Q  has  a  canonical  symplectic 
structure  w  =  ~d6  where  6  is  an  intrinsically  defined  one-form,  i.e.,  in  the  coordinate 
notation  &  —  £,•  p,dg,-. 

2.2.2  Momentum  Mapping 

Consider  two  symplectic  manifolds  (Af,w)  and  ( N,p ).  A  map  /  :  ( M,u )  ( N,p ) 
is  symplectic  if  /*  •  p  =  u>.  Let  (P,w)  be  a  connected  symplectic  manifold.  An 
action  9  :  G  X  P  P  of  a  group  G  is  a  symplectic  action  if  for  each  g  6  G  the  map 
9t  :  P  — ►  P  :  x  i— i ►  is  symplectic.  Let  g  and  g*  be  the  Lie  algebra  of  G  and 

its  dual  respectively.  A  mapping  J  :  P  —*  g*  is  a  momentum  mapping  of  the  action 
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$e  if  for  every  £  €  g 


«*P(0)(*)  = 

where  J(^)(x)  =  J(x)  •  £  and  is  and  inner  product  s. 

Let  $  be  an  action  on  Q  and  let  $T*  =  T*$g- 1  be  the  lifted  action  on  P  =  T*Q. 
Then  this  action  is  symplectic  with  respect  to  the  canonical  symplectic  structure  and 
has  an  Adm  equivariant  momentum  mapping  given  by 

J  :  P  —*  g*  :  J(p)  =  J(p).£  =  p*£c(g) 

where  £  E  g  and  (q  is  the  infinitesimal  generator  of  the  action  corresponding  to  £. 
The  momentum  for  a  vector  field  X  on  Q  is 

P{X)  :  T*Q  ->  R  :  p  p  •  X(g) 

The  following  diagram  commutes. 

.  0*  w  . 

T  Q - ►  T  Q 


If  L  is  the  Lagrangian  on  TQ  with  =  ( FL)*9  and  L  is  invariant  under  the 
action  =  T$qi  then  the  momentum  mapping 

?(£)(«,)  =  FL(v,)-(oW 

is  Ad*  equivariant  and  J  is  an  integral  of  the  Lagrangian  equations. 

EXAMPLE :  Consider  Q  =  Rs  and  G  =  50(3);  one  can  identify  so(3)  =  TtSO{ 3), 
the  Lie  Algebra  of  50(3),  with  Rs  as  follows.  Define 

0  —  Xi  xj 

Rs->rt50(3)  :  x  =  [xx,x2,x3]T  x  =  xs  0  -xj 

— x2  xj  0 

3see  Abraham  and  Marsden  [2]  page  115  for  definition  of  the  inner  product 
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An  easy  computation  shows  that  for  x,y  £  Rs 

(xxy)  =  [x,y]  =  x-y-y-x. 

Thus  TeSO( 3)  =  so(3)  the  Lie  algebra  of  50(3)  can  be  viewed  as  Rs  with  the  cross 
product  being  identified  with  the  Lie  bracket. 

Consider  the  group  action  $  :  50(3)  X  Q  — ►  Q  :  {A,q)  i->  Aq.  If  £  £  Rs  then 

A 

£  €  so(3).  Now  the  infinitesimal  generator  £q  can  be  seen  to  be 

fo  =  Jt  (*(exp(<?),  ?))  |t=0 

=  jt  (exp («?)  •  ?)  |(=0 

=  £  ■  9 
=  £  X  q 

The  momentum  mapping  J(p)  where  p  £  T*Q  is  given  by 

=  P'tq 

Working  on  TQ  we  have 

^(?)(?,v)  =  <v,iq> 

=  <  q  x  v,  £  > 


or 


J(?,u)  =  q  xv 
which  is  just  the  angular  momentum  of  the  system. 
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CHAPTER 

_ THREE 

PLANAR  MULTIBODY  SYSTEMS 


The  problem  of  dynamics  and  control  of  multibody  systems  has  long  been  of  spe¬ 
cial  interest  to  the  aerospace  community  [12]  [22]  [47]  [57].  The  methods  used  for  the 
generation  of  equations  of  motion,  in  general,  could  be  grouped  broadly  under  vec¬ 
torial  mechanics  (Newton-Euler  formulation)  and  analytical  mechanics  (Lagrangian- 
Hamiltonian  formulation)  [38]. 

The  Newton-Euler  method  involves  the  application  of  Newton’s  laws  of  motion  and 
Euler’s  equations  to  generate  the  necessary  dynamical  equations.  Since  the  method 
utilizes  physically  observable  quantities,  for  relatively  small  systems,  it  is  intuitively 
appealing.  However,  for  more  complex  systems  the  computations  can  beGome  tedious. 
Efforts  have  been  made  to  take  advantage  of  ways  to  mechanize  the  formulations  [22], 
[13],  [59].  Symbolic  processing  capabilities  have  been  exploited  to  generate  automat¬ 
ically  the  dynamical  equations  [14],  [49],  [60]  .  More  on  this  in  Chapter  7. 

The  Lagrangian  method  requires  the  formulation  of  a  scalar  quantity  called  the 
Lagrangian.  This  quantity  is  the  difference  between  the  kinetic  and  potential  ener¬ 
gies  of  the  system.  A  systematic  procedure  of  partial  differentiation  involving  the 
Lagrangian  would  result  in  a  set  of  second-order  ordinary  differential  equations  which 
represent  the  dynamics  of  the  system  [38]. 

Hamiltonian  methods  involve  the  formulation  of  a  single  variational  equation  using 
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Hamilton’s  principle  [l].  Much  of  Newtonian  mechanics  is  inherent  in  this  single 
equation.  The  representation  of  the  dynamical  system  is  by  a  set  of  canonical,  first 
order,  ordinary  differential  equations  with  a  simple  structure.  An  underlying  Poisson 
structure  provides  an  elegant  representation  of  the  first  order  dynamical  equations  in 
terms  of  the  bracket.  The  Poisson  structure  can  also  help  us  identify  the  kinematically 
conserved  quantities,  the  Casimirs.  Casimirs  could  be  effectively  employed  to  study 
the  stability  of  equilibria.  The  chaotic  solutions  produced  by  the  splitting  of  the 
homoclinic  orbits  when  the  system  is  perturbed,  could  be  studied  using  Melnikov’s 
method  (see  Guckenheimer  and  Holmes  [17]).  Hamiltonian  methods  have  been  applied 
to  dynamical  systems  in  space  in  the  past  by  Pringle  [44],  [45]  and  more  recently 
by  Krishnaprasad  [32],  Krishnaprasad  and  Berenstein  [33],  and  Krishnaprasad  and 
Marsden  [34]. 

Satellite  Dynamicists  are  interested  in  both  rigid  and  flexible  body  formulations 
[38].  In  this  research  we  will  focus  on  satellites  modeled  as  rigid  bodies. 

Here  we  formulate  the  dynamics  of  a  planar  system  of  arbitrary  number  of  rigid 
bodies  connected  in  the  form  of  a  tree  structure  using  Lagrangian  and  Hamiltonian 
formulations.  We  start  by  giving  the  notation  used,  followed  by  the  Lagrangian 
dynamics.  The  Legendre  transformation  which  maps  the  Lagrangian  equations  to 
the  Hamiltonian  equations  is  given  next.  We  discuss  the  underlying  symmetries  and 
the  procedure  of  reduction  and  effectively  use  it  to  get  the  reduced  dynamics  in  the 
Hamiltonian  setting. 

3.1  Notation 

Consider  a  system  of  arbitrary  number  of  planar  rigid  bodies  connected  in  the 
form  of  a  topological  tree1  (no  closed  paths).  To  define  the  system  mathematically, 
each  individual  body  should  be  assigned  a  unique  label.  A  set  of  consecutive  integers 
is  used  as  labels.  Let, 

N  -  total  number  of  bodies  in  the  system. 

1i.e.,  there  exists  a  unique  path  from  one  body  to  another  -  also  referred  as  an  open  kinematic 
chain 
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Body  1  could  be  chosen  arbitrarily.  Let  body  1  be  the  principal  or  the  main  body 
of  the  system.  All  other  bodies  are  given  distinct  integer  labels  varying  from  2  to 
N  inclusive.  For  simplicity  of  computation,  it  is  convenient  to  label  the  bodies  such 
that  the  body  labels  are  of  increasing  magnitude  along  any  topological  path  starting 
at  body  1. 

A  body  is  considered  inboard  of  body  »  if  it  is  contained  in  the  topological  path 
connecting  body  1  and  body  «  .  A  body  is  considered  outboard  of  body  »  if  it  is  not 
inboard  of  body  *  ,  and  if  body  i  is  contained  in  the  topological  path  connecting  body 
1  and  the  body  under  consideration. 

Similarly  it  is  necessary  to  assign  unique  labels  to  all  joints.  An  efficient  way  of 
doing  it  is  the  following.  The  joint  connecting  body  i  and  the  body  which  is  contiguous 
to  and  inboard  to  body  i  is  labeled  i-1.  Hinge  0  is  assumed  to  be  the  center  of  mass 
of  body  1. 

Let, 

i  -  body  label  of  body  *  , 

J{i)  -  body  label  of  the  body  contiguous  to  and  inboard  of  body  i  , 

i  -  joint  label  of  joint  connecting  body  »’  and  body  J(t);  also  referred  to 

as  joint  (t  —  1), 

Sk-i,i  -  set  of  all  body  labels  associated  with  those  bodies  lying  on  the  topo¬ 
logical  path  from  joint  k  —  1  to  the  center  of  mass  of  body  t  , 

Figure  3.1  shows  one  way  of  labeling  a  system  of  12  bodies  connected  in  the  form 
of  a  tree.  The  topological  path  connecting  joint  0  (body  1)  and  body  6  ,  i.e.,  the  set 
of  body  labels  in  the  path  connecting  body  1  (joint  0)  and  body  6,  50,6»  is  given  by 

S0l6  =  [1,2,  5, 6].  (1.1) 

Here,  body  5  is  considered  inboard  of  body  6.  Body  6  is  considered  outboard  of  bodies 
1,  2  and  5.  Body  9  is  neither  inboard  nor  outboard  of  body  5.  The  body  contiguous 
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Figure  3.1:  Planar  12-body  system  connected  in  the  form  of  a  tree 


to  and  inboard  of  body  5,  i.e.,  J(5)  is  body  2.  The  joint  connecting  body  5  and  J(5), 
(body  2)  is  labeled  joint  4  (i.e.,  joint  (5-1)). 

To  formulate  the  equations  of  motion  of  the  system  under  consideration,  define 
the  origin  of  the  inertial  coordinate  system  to  be  fixed  at  the  point  of  reference  with 
its  Z  axis  being  perpendicular  to  the  plane  of  the  paper.  Define  a  local  frame  of 
reference  for  each  body,  the  origin  of  which  is  fixed  to  the  body  center  of  mass.  The 
Z  axes  of  the  local  frames  of  reference  are  also  perpendicular  to  the  plane  of  paper. 
Refer  to  Figure  3.2. 

Oi  •  origin  of  the  local  frame  of  reference  of  body  t  located  at  the  center  of 
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mass  of  body  *  , 


00  -  center  of  mass  of  the  system  of  N  rigid  bodies, 

Oj  -  reference  point  and  origin  of  inertial  coordinate  system, 

Oi  -  angle  made  by  the  local  frame  of  reference  to  the  inertial  frame  of 

reference, 

Oij  •  angle  made  by  the  local  frame  of  reference  of  body  »  with  respect  to 
the  local  frame  of  reference  of  body  j\  equal  to  (0,-  —  0j), 

u>i  -  inertial  angular  velocity  of  body  *  along  the  Z  axis  (perpendicular  to 
the  plane  of  the  paper);  equal  to 

R(0i)  -  (2  x  2)  rotation  matrix  associated  with  body  i  , 

cos  (Oi)  —  sin(0,) 
sin(0,)  cos($,-) 

To  complete  the  description  of  the  multibody  system,  define  a  set  of  position 
vectors  (see  Figure  3.3)  and  associated  physical  parameters  of  the  system  :  let, 

aj  -  inertial  vector  from  joint  (*  —  1)  to  the  center  of  mass  of  body  i  , 

Pi  -  inertial  vector  from  joint  (J(i)  —  l)  to  joint  (i  —  1)  , 

on  -  vector  from  joint  (t  —  l)  to  the  center  of  mass  of  body  »  in  the  local 

frame  of  reference  of  body  *  , 

Pi  -  vector  from  joint  ( J (t)  —  1)  to  joint  (*  —  1)  in  the  local  frame  of  reference 
of  body  J(t), 

r  -  vector  from  the  reference  point  to  the  system  center  of  mass, 

T°  -  vector  from  the  center  of  mass  of  the  system  of  N  bodies  to  the  center 

of  mass  of  body  »  ,  ' 


(1.2) 
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body  j 

k 


Figure  3.3:  Body  and  inertial  vectors 

r<  -  vector  from  the  reference  point  to  the  center  of  mass  of  body  i  . 

See  Figure  3.2  for  a  description  of  some  of  these  vectors. 

Note  that  for  all  i  = 

r.  =  r  +  r°i 

=  r  +  rj  +  52  Pi  + 

/€S0,<.  J'/I 

—  +  52  Pi  +  ®»»  (1*3) 

jes0,,-.  y/i 

m,  -  mass  of  body  »  , 

N 

m  -  (!-4) 

IS  I 

-  moment  of  inertia  of  body  »  at  its  center  of  mass  and  along  an  axis 
perpendicular  to  the  plane  of  the  paper, 

Ki  -  kinetic  energy  of  body  i 


16 


3.2  Basic  Kinematics 


The  configuration  space  Q  for  the  planar  N-body  system  is  the  subset  of  Si? (2)  x 
•  •  •  x  SE{2)  ( N  copies  of  the  special  Euclidean  group  of  the  plane)  consisting  of  pairs 
((iZ(0i),ri), . . . ,  (i2(0N)»rAr))  and  satisfying  the  joint  constraints  : 


rk  =  ri  +  *52  Pi  + 

ieSo'k,  & i 

=  *i  +  R(0JU))0i  +  k  =  2 (2.1) 

ieSo't,  i*i 

Notice  that  Q  is  of  dimension  N  +  2  and  is  parametrized  by  0i, . . . ,  On  and,  say  n; 

i.e.  Q  «  S1  X  -  “S1  xlR2.  We  form  the  velocity  phase  space  TQ  and  momentum 
N  times 

phase  space  T*Q. 


The  Lagrangian  on  TQ  is  the  map  L  :  TQ  — ►  H,  obtained  by  summing  the  kinetic 
energies  of  the  individual  bodies,  i.e., 

L  =  (2.2) 

k-i 

where  if*  is  the  kinetic  energy  of  body  k.  The  equations  of  motion  then  are  the 

Euler-Lagrange  equations  for  this  L  on  TQ  : 

’  d_  fdL\  _  dL_ 
dt  \dg)  dq 

where,  q  =  [^x, . . . ,  fl/v]7-  Equivalently,  they  are  Hamilton’s  equations  for  the 

corresponding  Hamiltonian. 

To  formulate  L  in  terms  of  the  kinematic  parameters  we  proceed  as  follows  :  let 
X*  denote  a  position  vector  in  body  k  relative  to  the  center  of  mass  of  body  k,  and  let 
Pk(Xk)  denote  the  mass  density  of  body  k  (see  Figure  3.4).  Then  the  inertial  position 
of  the  point  with  material  label  X*  is 


x*  =  R  (0*)X*  +  r*. 


Thus 


xk  =  R{Ok)yLk  +  r*, 
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body  k 


Figure  3.4:  Position  vectors  in  body  t 
and  so  the  kinetic  energy  of  body  k  is 

Kt  =  \  11  **  n!  *** 

=  i  /  />»(X»)  <  R(h)y.k  +  *»,  A(**)X»  +  h>  d'xt 

Z  JBk 

-  J  f  Pk(Xk)[<  R{0k)XktR{0k)Xk  >  +  2  <  R(0k)Xk,tk  > 

Z  JBk 

+  II  K  IIVX*.  (2.4) 

Here  Bk  denotes  the  extent  of  body  k  and  pkd?Xk  denotes  its  mass  measure. 

But 

<R{0k)Xk,R{6k)Xk>  =  tr(R(0k)Xk(R(0k)Xk)T) 

=  tr(R(0k)XkXZR(0k)T),  (2.5) 

and 

f  Pk(Xk)  <  R(0k)Xkirk  >  d*Xk  =  <  R(0k)  f  pk(Xk)Xkd3Xkyrk  >  =  0, 

J$k  JBk 

(2.6) 

since  X*  is  the  vector  relative  to  the  center  of  mass  of  body  k.  Substituting  (2.5)  and 
(2.6)  into  (2.4)  and  defining  the  matrix 

I*  =  f  p(Xk)XkXld3Xk,  (2.7) 
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we  get 


K„  =  i|r(A(*t)l‘tf  («*)*•)  +  lmt  ||  i „  ||* 

=  itr{R(h)liR(ei)T)  +  imi  ||  rj  +  r  ||* 

=  j<r(fl(».)I*A(».)r)  +  im4[||rj||*  +  2<r;,r>  +  ||  r  ||>] .  (2.8) 


For  later  convenience,  we  shall  rewrite  the  energy  (2.2)  in  terms  of  u>k  =  and 
rj’s.  To  do  this,  note  by  definition  the  vector  from  the  reference  point  to  the  center 
of  mass  of  the  system  in  the  inertial  coordinate  system  is  given  by 


1  N 

r  =  ~  E  m*r*- 

m  . 


k  =  1 


Substituting  for  r,  from  (1.3)  and  simplifying  we  have, 

1  N 

r  =  ri  +  —  £  mt 
m 


1=2 


y  &}  + 

jeso,,.  j*  i 


Keeping  in  mind  that  rx  =  r  +  rj,  we  rewrite  the  above  equation  as, 


1  N 

ri  =  -  —  e 

m  (=2 


XT  Pi  +  ai 

Ue50,{>  ;/l 


where  €j  = 


N 

1=2 
mi 


y  Pj  + 

L;e50,(,  ;Vi 


m 


We  know  that, 


(2.9) 


(2.10) 


(2.11) 

(2.12) 


rJb  —  rX  +  E  Pi  +  ak • 

yeJo,*.  j'/i 


From  (2.11)  and  (2.13)  we  get, 


N 


r°k  =  ~E  c< 


y  Pt  + 


1=2  lies o,(l  .yi 


+  E  Pi  + 

ye-50,t, 


* 


—  E  (a*.iA  +  , 

<  =  i 


(2.13) 


(2.14) 
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where, 


hj  = 


Ki  = 


1  -  EyL2  jQ  1  S  S’o.t  and  fc  7^  1, 
(  -  EyLj  iij«l  otherwise  , 

1  if  leSoj , 

0  otherwise, 

1  —  cj  if  l  =  k, 

-  £,  if  /  #  Jfc. 

Differentiating  (2.14)  w.r.t.  time  we  get 

AT 


r't  —  £  a*,iA  +  • 


(  =  1 


Now  we  substitute 


into  (2.14)  to  get 


Pi  =  R{0j(i))$i  so  Pi  =  R{9j[i))Pi , 

aj  =  JE(0j)3|  so  dj  =  R(0i)ai, 


N 

ri  =  £  \ak,iR{9j[i))Pi  +  &jfc,j.R(0i)aj] 

= 

1=1 


and  into  (2.15)  to  get 


N 

K  =  £  [ak,,R(0j(i))Pi  +  &M^(^)Si] 

la  1 

1=1 


where  £*,1  is  given  by 


(2.15) 


(2.16) 


(2.17) 


$k,i  =  £  [a*, ,•/?<]  +  bk, i  an 

Vf,  *.«.  J(i)  =  1 


%  J 


(2.18) 
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Also  from  (1.2) 


R(0i)  = 


d_ 

dt 


cos  (0i)  —  sin(0/) 
sin(^j)  cos(0|) 

—  sin(^)  —  cos(^) 
cos(0j)  —  sin(0j)  J 


ui  :=  R{0i ) 


0  —u>i 
u>i  0 


:=  R{6i)Qi. 

From  (2.17)  and  (2.19)  we  find 

ii r**  ii 2  =  <  uxtofoji  n*(*i)2U> 

y=i  i=i 

j=l 1=1 

=  Eiia&nx 

y=i 


N-l  N 

+  *E  E 

y=i  i=y+i 


cos  (0y  —  0i )  sin(0y  —  $i) 

—  sin(0y  —  0j)  cos  {0  j  —  0j) 


(2.19) 


.4(/cj,o;_,^2.20) 


Substituting  from  (2.18)  for  6k, i  in  the  above  equation  we  get, 

«f:ir  =  eii^ii’"‘  +  2E  £  (piA  +  cost#,  - «,) 

y=x  y=i  J=y+i 


+  p IJh  -  s'm(ei  -  0')}  <*>iwy 


=  E II  II2  ui  +  2  E  E  {p*A]  cos(«,-,) 

y=i  y=x  J=y+i 


+  [l&j  x  hA}  sin(^y,i)}  UiUj, 


(2.21) 


where  6j,i  =  (0;  —  $i )  and 


l^.y  x  ^u!  :=  ptj^i,i —  • 


(2.22) 


Recall  that  the  Lagrangian  for  the  system  which  is  also  the  total  kinetic  energy  of  the 
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system  is  then  given  by, 

L  =  |  E  tr{R{et)lkR[Bt)T)  +  |  fm,  [||  rj  ||s  +  2  <  rj.r  >  +  ||  r  f]  .  (2.23) 


i=l 


But 


k=l 

JV 

X]  mkr°k  =  0  so  =  °> 

k=l  k=  1 


(2.24) 


and  from  (2.19)  we  can  see  that 

tr(R(dk).lk.R(0k)T)  =  tr(R{dk)T.ullk.Qk.R{9k))  =  tr{C%  Ifcu>k) 

=  tr(ul  u>kIk)  =  tr(  *  I*)  =  u^tr(I*) 

0  u* 

■=  <±>lh>  (2.25) 

where, 

h  =  f  p{Xk,Yk){X\  +  Y?)dXkdYk, 

J  &h 

is  the  moment  of  inertia  of  the  body  k  about  the  Z  axis  of  the  local  frame  of  reference. 
Substituting  (2.24)  and  (2.25)  into  (2.23)  gives 

L  =  \  Zj  {7*wk  +  ™k  ||  r°k  ||2  +  mk  ||  r  ||2}  . 

Finally  simplifying  using  (2.21)  we  get 

L  = 


=  +  *E  E  XyiWu)"*".}  +  im||rp 

z  (k=i  ;=i  1=3+1  )  1 


-  Vj  U3  +  11  P  11 

-  2-  J-  +  ~2m~* 


where,  p  =  mr  is  the  system  linear  momentum ,  m  is  as  given  by  (1.4), 

W.  —  [wi,***,a;jv]r , 


(2.26) 


(2.27) 


J  = 


h 

A12 


Ai2 

h 


A IS 


'2  N 


n  %2N  •••  In 


(2.28) 
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is  a  symmetric  nonsingular  matrix;  for  the  sake  of  clarity  the  arguments  of  A ji  have 
been  dropped, 


mjk  (%kj-6kl) 


cos(0u) 


+ 


[]C  mk(\6k:  X  £*i|) 


U=1 


sin(0u),(2.29) 


note  that  A  ji  =  A /,•; 


Ik 


Ik 


+ 


Z  mi  II  II2  • 


3=1 


(2.30) 


<w 

Remark  3.2.1  :  The  Ik  of  a  body  k  is  the  moment  of  inertia  of  an  “augmented” 
body  k  with  lumped  reduced  masses  at  the  joints  through  which  it  is  connected  to 
other  bodies.  It  is  to  be  noted  here  that  the  concept  of  “augmented”  bodies  arises 
naturally  in  the  context  of  multibody  problems  in  our  formulation. 


Remark  3.2.2  :  The  term  “augmented  body”  was  first  used  by  Hooker  and 

Margulies  [22]  in  the  context  of  tree  connected  multibody  dynamics.  They  define 
“augmented  body”  k  as  follows.  Consider  a  body  k  in  the  multibody  system.  For 
each  joint  (/  —  1)  on  the  body  k  the  mass  distribution  of  body  k  is  augmented  by  the 
mass  equal  to 

mk  =  Z  mr 

3  €St_,,,t  jjik 

Obviously  the  mass  of  the  “augmented  body”  is  defined  as  equal  to 

mk  =  mk  +  mk- 

Vi.  i.t.  J{l)=k 

The  center  of  mass  of  the  “augmented  body”  (also  known  as  the  connection  barycen- 
ter  or  simply  barycenter)  could  be  obtained  easily.  The  inertia  dyadic  for  this  “aug¬ 
mented  body”  about  its  barycenter  figures  prominently  in  the  formulation  of  the 
system  dynamics.  Also  see  Wittenburg  [59].  In  contrast  a  quick  look  at  (2.30)  and 
(2.18)  indicates  that  in  our  formulation  we  use  lumped  reduced  masses  at  the  joints 
instead  of  lumped  masses. 
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3.3  Legendre  Transformation 


The  Legendre  Transformation  maps  the  Lagrangian  equations  into  the  corre¬ 
sponding  Hamiltonian  equations.  Classically,  the  name  Legendre  transformation  (see 
Courant  and  Hilbert  [8])  is  reserved  for  a  map  from 

£(2,£)  =  £r£  ~  £(£>2)i  (3-l) 

where  L  :  TQ  — *  R  is  the  Lagrangian,  H  :  T*Q  -+Ris  the  Hamiltonian,  and 

£  =  [9i . 9w]r,  (3.2) 

t  ,t  dL 

P  =  [Pi,...,Ptf]  =  -jp-.  (3.3) 

Here,  we  state  a  theorem  which  gives  the  Legendre  transformation  explicitly  with 
the  proof  given  in  coordinate  notation.  Using  this  transformation  we  construct  the 
Hamiltonian  for  the  system  of  N  planar  rigid  bodies  connected  in  the  form  of  tree 
structure.  Some  definitions  are  in  order  before  we  give  the  theorem. 

Definition  3.3.1  :  Let  ( E,irtM )  and  ( Ftp,M )  be  fiber  bundles  with  the  same 
base  space  M.  Let  it  :  E  — »  M  :  (nm  and  p  :  F  —►  M.  If  e0  =  7r_1m  C  E , 
then  Cq  is  a  fiber  of  E.  The  concept  of  fiber  bundle  generalizes  the  notion  of  a  product 
space. 

Let  /  :  E  —*  F  be  a  smooth  map  which  is  fiber  preserving,  then  ?r(e)  =  po  /(e), 
i.e.,  the  following  diagram  commutes. 


-►  F 


Let  Jr(em)  =  m  £  M  and  fm  =  f\Em  =  /|7r-1(m).  Then  the  fiber  derivative 
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of  fm  (Frechet  derivative  of  fm)  evaluated  at  em  is  the  map 

Ff :  E  — ►  U  L{EmtFm), 

m  6  Si 

Cm  »-►  Dfm(em)  €  L(BmtFm). 


For  our  purposes  we  specialize  M  to  be  the  configuration  space  Q,  E  to  be  TQ 
the  tangent  bundle  of  Q,  and  F  to  be  Q  X  R.  We  consider  the  Lagrangian  L  on  the 
velocity  phase  space  TQ  as  the  map  L  :  TQ  —>  1 1.  Let  L  :TQ  — ►  QxR  :  t/9  gxR. 

AT 

Then  the  ^6cr  derivative  FL  is 

FL{vq)  «  FL(u9)  =  £>L9  e  L(T(??,R)  =  T’Q 

The  Hamiltonian  on  the  momentum  phase  space  T*Q  is  a  map  H  :  T*Q  —*■  R  that 
can  be  expressed  in  terms  of  FL. 

Definition  3.3.2  :  The  action  A  is  defined  as  A  :  TQ  — ►  R,  by  A{vz)  =  FL(yx).vx\ 
in  coordinate  notation  A(q,q)  —  £Tff-  The  energy  function Et  defined  by  E  —  A  — L 
is  determined  on  TQ  by  L. 

When  E  is  translated  to  T*Q  by  means  of  the  fiber  derivative  of  L,  FL  :  TQ  — ► 
T*Q,  we  get  a  suitable  Hamiltonian.  Then  solution  curves  in  both  T*Q  and  TQ  will 
coincide  when  projected  to  Q,  i.eM  the  solution  curves  of  the  Hamiltonian  and  the 
Lagrangian  equations  coincide  on  Q.  Tq  and  Tq  are  maps  from  TQ  and  T‘Q  to  Q 
respectively.  The  following  diagram  helps  keep  the  relations  straight. 

R 

T*Q  "« - TQ - - - ►R 

T 

Q 


Definition  3.3.3  :  Let  M  be  a  manifold.  Define  tm  :  TM  — +  M  and  c(f)  :  I  — ♦ 
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TM.  Then  tm  o  c  :  I  — *  M  is  a  base  integral  curve. 

Theorem  3.3.1 :  Let  L  be  the  Lagrangian  on  Q  with  a  diffeomorphic  fiber  derivative 
FL  :  TQ  — ►  T*Q ,  and  let  H  =  E  o  ( FL)~X  :  T*Q  —*  R,  where  E  is  the  energy  of 
L.  Then  the  vector  fields  XE  and  XH  are  related  by  {FL),XE  =  Xu.  Furthermore 
XE  and  Xh  have  the  same  base  integral  curves. 

Proof  :  The  proof  of  this  theorem  is  given  here  in  coordinate  notation  (for  more 
details  and  a  ‘coordinate  free’  proof  see  Theorem  3.6.2  and  Example  3.6.10,  Abraham 
and  Marsden  [2]). 

Let  q  =  [?x*  •  •  • » ?Ar]r  and  p  =  [pi,... ,  p^]r.  By  the  definition  of  fiber  derivative 
FL  we  have 


FL-.TQ  -  T'Q  :  (j.j)  «  (,,£)  =  (j,  ||  V 

Recall  from  the  definition  of  the  action  we  have  A(q,q)  =  qTjr  and  p  =  ||  from 
(3.3),  so, 


=  Eo(FL)  *(£,£) 

=  ^(£,£)  ~  L(q,q) 
dL 

=  iaj  -  Mi- s) 

=  qTp  -  L(q,q),  (3.4) 


which  is  the  Hamiltonian. 
Taking  partials  on  (3.4)  yields 


dH 

dpi 


=  ft  + 


-  £ 


dL  dq, 
dq,  dpj 


ft- 


Keeping  in  mind  the  Lagrangian  dynamical  equations 
we  find  that 


(3.5) 

If  =  0.  »nd  (3.3) 
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or  from  (3.4)  we  have, 


dH 

dq{ 


_dL 

dqi 

~Pi- 


(3.6) 


So  from  (3.5)  and  (3.6)  we  find  that 


•  _  Off 

2  ~  dp 

i  _  BH 

E  -  a,  ) 


(3.7) 


which  are  the  Hamiltonian  dynamical  equations  with  the  Hamiltonian  given  by  (3.4). 
Q.E.D. 


3.3.1  Hamiltonian  for  the  Multibody  Problem 


In  the  planar  multibody  problem,  the  configuration  space  Q  can  be  coordinatized 
by  (0i,...,0tf,ri)  and  the  tangent  space  TQ  by  {0U . . . ,  0N,  u>i, . . .  ,vN),  i.e.,  £  = 
[wj, . . .  ,cj//]r.  The  Lagrangian,  from  (2.26),  is  given  by 


1  T  II  P  II2 

L  —  —oj  Ju>  +  — — - . 

2~  2m 

With  p  =  /z  =  [ii\, . . .  ,Miv]T>  a  simple  calculation  from  (3.3)  indicates  that 

dL 


(3.8) 


(*  = 


du 
=  J  v. 


(3.9) 


Constructing  the  Hamiltonian  analogous  to  (3.1)  using  (3.9)  leads  to 

B  =  |/£TJ'V  +  (3.10) 

where  n  =  [/*i , . . .  ,A*Ar]r  and  (jl^s ,  k  =  1,  •  •  •  ,JV,  axe  the  conjugate  momenta. 


Remark  3.3.1:  By  making  use  of  the  reduction  technique  due  to  Arnold  [l]  and 
Marsden  and  Weinstein  [41]  we  can  prove  that  in  the  general  case  of  a  system  of  N 
planar  rigid  bodies  connected  in  the  form  of  a  tree,  the  Hamiltonian  is  indeed  given 
by  (3.10),  with  the  linear  momentum  equal  to  zero.  This  we  do  in  the  next  section. 
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3.4  Symmetries 

In  this  section  we  recognize  the  basic  symmetries  involved  in  the  planar  N  rigid 
body  problem  and  reduce  the  dynamics  accordingly.  The  reduction  technique  used 
here  goes  back  to  Arnold  [l]  and  Marsden  and  Weinstein  [41].  One  starts  with 
a  Poisson  manifold  P  and  a  Lie  group  G  acting  on  P  by  canonical  transforma¬ 
tions.  The  reduced  phase  space  P/G  (assume  it  has  no  singularities)  has  a  natural 
Poisson  structure  whose  symplectic  leaves  are  the  Marsden-  Weinstein-Meyer  spaces 
J~1(fi)/Gtt  «  J~l(0)/G  where  p,  €  g%  the  dual  of  the  Lie  algebra  of  G,  J  :  P  — ►  g* 
is  an  equivariant  momentum  map  for  the  action  of  G  on  P,  G^  is  the  isotropy  group  of 
p  (relative  to  the  coadjoint  action)  i.e.,  (7M  =  {0  6  G  :  Ad*-lp  =  /zj,  and  0  is  the 
coadjoint  orbit  through  p.  The  coadjoint  orbit  0,  is  even  dimensional  and  possess  a 
natural  symplectic  structure  known  as  Kirillov  2-form  (see  Kirillov  [35]).  For  details, 
we  refer  the  reader  to  Chapter  2  for  a  brief  exposition  on  momentum  map  (also  see 
Abraham  and  Marsden  [2]  Chapter  4,  Section  4.3)  and  Marsden  and  Weinstein  [41] 
for  the  Marsden- Weinstein  reduction  theorem. 

3.4.1  Reduction  to  the  Center  of  Mass  Frame 

We  reduce  the  dynamics  by  the  action  of  the  translation  group  1R.2.  This  group 
acts  on  the  original  configuration  space  Q  by 

v  •  ((R(9i),t1),...,(R(0n),tn))  =  ((P(0i),ri  +  \),...,(R(0N)frN  +  v))  (4.1) 

This  is  well-defined  since  the  joint  constraints  (2.1)  are  preserved  by  this  action.  The 
induced  momentum  map  on  TQ  is  calculated  by  the  standard  formula 

A  =  (4.2) 

or  on  T*Q  by 

=  P»  £<?(?)>  (4-3) 

where  £*  is  the  infinitesimal  generator  of  the  action  on  Q.  (see  Abraham  and  Marsden 
[2]).  To  compute  (4.2)  and  (4.3)  we  parametrize  Q  by  9i, . . . ,  9 ^  and  r  with  r*5  k  = 
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1 ,N  determined  by  (2.14).  From  (2.26)  we  see  that  the  momentum  conjugate  to 

r  is 

dL 


and  so  (4.2)  and  (4.3)  gives 


p  =  di  =  mr’ 


J(  =  <P,£>, 


(4.4) 


Thus  J  =  p  is  conserved  since  H  is  cyclic  in  r  and  so  is  translation  invariant.  The 
corresponding  reduced  space  is  obtained  by  fixing  p  =  Po  and  letting 


-Pp.  =  J-'(Po)/E!. 


(4.5) 


(see  Chapter  4,  Abraham  and  Marsden,  [2],).  But  PPo  is  clearly  isomorphic  to 
T'fS1  X  •••  X  S1)  i.e.  to  the  space  of  0i,...t  Os  and  their  conjugate  momenta 

N  timet 

(//i, . . .  ,Hs)'  The  reduced  Hamiltonian  is  simply  the  Hamiltonian  corresponding 
to  (2.26)  with  p  regarded  as  a  constant. 

Note  that  in  this  case  the  reduced  symplectic  manifold  is  a  cotangent  bundle,  in 
agreement  with  the  cotangent  bundle  reduction  theorem  (Abraham  and  Marsden  [2], 
Kummer  [36]).  The  reduced  phase  space  has  the  canonical  symplectic  form;  one  can 
also  check  this  directly  here. 

In  (2.26)  we  can  adjust  L  by  a  constant  and  thus  assume  p  =  0;  this  obviously 
does  not  affect  the  equations  of  motion. 

Reduced  Hamiltonian 


The  Lagrangian  from  (2.26)  is  given  by 

1 


Q  2 

L  —  —  Jw  +  — — - . 

2  ~  2m 


(4.6) 


Assuming  p  =  0,  and  using  the  Legendre  transform  and  the  corresponding  conjugate 
momentum  vector  (see  section  3.3.1)  j±  =  3 u,  we  get  the  Hamiltonian  as, 


H  = 

2-  - 


(4.7) 
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It  is  to  be  noted  here  that  J  =  J(g)  where  £  =  [0i,  •  •  • ,  0^  ]. 

Let  us  observe  that  the  reduced  system  is  given  by  geodesic  flow  on  S1  X •  •  •  x  51 

- - - — 

N  timet 

since  (2.26)  is  quadratic  in  the  velocities.  Indeed  the  metric  tensor  is  just  the  matrix 
J  given  by  (2.28). 

Remark  3.4.1  :  The  reduction  to  center  of  mass  coordinates  here  is  somewhat 
simpler  and  more  symmetric  than  the  Jacobi-Haretu  reduction  to  center  of  mass 
coordinates  for  n  point  masses.  (Just  taking  the  positions  relative  to  the  center  of 
mass  does  not  achieve  this  since  this  does  not  reduce  the  dimension  at  all!)  What  is 
different  here  is  that  the  bodies  are  connected  by  revolute  joints,  and  so  by  (2.13), 
rj^’s  for  j  =  1, . . . ,  N,  are  determined  by  other  data. 

3.4.2  Reduction  by  Rotations 

Recall  from  the  previous  section  that  the  configuration  space  after  reduction  to 
the  center  of  mass  frame  is 

Q  =  [S1  x  ...  x  S1), 

N  timet 

and  so  the  cotangent  space  is  T *(5*  X  »»•  X  Sl). 

N  timet 

There  exists  a  rotational  symmetry  group  5 1  which  acts  on  the  cotangent  space 
as  follows  : 


0  •  ((^i» /*i)»  *  *  *  >  (^jv, Mat))  =  ((0i  +  0,nx),'  •  • ,  (0N  +  0,/ijv)) 

=  ((0ii^i)»'”»(0Ar.MJv))- 


To  complete  the  reduction,  we  reduce  by  the  diagonal  action  of  S 1  on  the  config¬ 
uration  space  S 1  X  •  •  •  X  S 1  that  was  obtained  in  Section  3.2.1.  The  momentum  map 

N  «  tmmt* 

S  timet 

for  this  action  is  found  next.  We  know  that  S1  is  diffeomorphic  to  SO (2)  (the  special 
orthogonal  group  of  (2  X  2)  matrices  i.e.,  A  6  50(2),  AAT  =  /,  and  det(vl)  =  1). 
The  Lie  algebra  of  50(2)  is  so(2)  (  skew-symmetric  matrices  with  determinant  equal 
to  k,  k£0),  i.e.,  TeSO( 2)  =  so( 2). 
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For 


£  = 


0  1 

-1  0 


where  £  6  so(2),  the  exponential  of  £  is 


exp(f£)  =  exp 


t 

\  L 


0  1 
-1  0 

cos(£)  sin(£) 

—  sin(f)  cos(t) 

So  exp(f£)  €  0(2). 

The  infinitesimal  generator  $q(q)  can  now  be  calculated  as  follows 
£<?(?)  =  ^(exp(tO»9)  It — o 


d.( 
=  dt*( 


cos  (£)  sin(t) 

—  sin(t)  cos(f) 

-  +t,...,0N  +t)  |(_0 

=  (1 . 1). 


,  {0U...,0N)  )  |t=0 


The  momentum  map  is  given  by 


J:T'Q  -  g*  =  so* (2), 
with  the  momentum  P  :  TQ  — ►  R 

PM  =  *(0(*tJ 

=  FL{vq)'  tQ(q) 

=  (vQ,tQ(q))  on  TQ. 

Note  that  the  metric  here  is  the  Reimannian  metric  and  the  inner  product  *{,)’  is 
given  by  (z,  y)  =  xT Jy. 

P{?q)  =  ((^l****  >uat)  >  (!>••■  .1)) 

=  [l,...,l]Ju/  on  TQ, 
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or  on  T’Q,  for  a,  6  T*Q  we  have 

P(aq)  =  Hi  +  --•  +  Hn- 

i.eM 

•J((0i»A*i)»***»(0.v»Mjv))  =  Mi  +•••+  Hn.  (4.8) 

For  purposes  of  later  stability  calculations,  we  shall  find  it  convenient  to  form  the 
Poisson  reduced  space 

P  :=  T'{Sl  x--‘XSl)fSl  (4.9) 

N  tim* » 

whose  symplectic  leaves  are  the  reduced  symplectic  manifolds 

Pn  =  J~1(y)/S1  C  P  (4.10) 

We  coordinatize  P  by  #*,/{*),  k  =  2, . . . ,  N,  and  y y,  j  =  1, . . .  JV;  topologically, 

P  =  S1  x  •••  x  Sl  xR*.  (4.11) 

(AT-l)  timtt 

The  Poisson  structure  on  P  is  computed  in  the  standard  way:  take  two  functions 
F{9itJWtk  ==  1, •  •  • , N, yiy . . . , i*n)  and  H{0k,J(k)tk  —  2 ,...,JV,  A^i,  < . .  ,/£w)<  Re¬ 
gard  them  as  functions  of  9\, . . .  ,0jv,A*i,  . . .  ,/ijv  by  substituting  Ok>j(k)  =  Ok  —  0/(4) 
and  compute  the  canonical  bracket. 


Lemma  3.4.1  :  The  Poisson  structure  on  P  is  given  by  the  non-canonical  bracket 


{/,<?}  = 

TV 

E 

4=2 

ff a/ 

df' 

\  dg 

f  dg 

dg  \ 

df  1 

,(4.12) 

dyk, 

/  ^^4,7(4) 

\dyj(k) 

dyk) 

dOk,j{k). 

where  /,  :  P  •— i ► 

R. 

Proof  :  The  canonical  bracket 

on  the  Poisson  manifold  P  is  given  by  (see 

Chapter 

2): 

{/.*>  = 

tf 

E 

4=1 

df  dg 

dg  df ' 

(4.13) 

dyk  dOk 

dyk  86 k 

* 

Since  f,g  :  P  H  and  P  is  coordinatized  by  yi, . . .  k  =  2 ,  ...,jV, 

we  have  /  and  being  functions  of  A*i>  •  •  • ,  (0*  —  0;(t)),  &  =  2,  ...,JV  ,  since 

0M(*)  =  (0*  -  0/(4))- 


32 


Therefore, 


dj_ 

dOk 

df 


df 


d9k,j{k)  ’ 

3/ 


(4.14) 

dOk,j{k)  ^  ^ 

Substituting  (4.14)  and  (4.15)  in  (4.13)  and  combining  suitable  terms  results  in  the 
non-canonical  bracket  (4.12).  Q.E.D.. 


CASIMIRS 


The  Casimirs  on  P  are  obtained  by  composing  J  with  Casimirs  on  the  dual  of  the 
Lie  algebra  of  S'1;  i.e.  with  arbitrary  functions  of  one  variable;  thus 

C  =  +  .  • .  +  Mjv)i  (4-16) 

results.  This  can  of  course  be  checked  directly. 

The  cotangent  bundle  reduction  theorem  asserts  in  this  case  that  the  reduction  of 
T"  (S1  x  •  •  •  x  S1)  by  S1  is  symplectically  diffeomorphic  to 

N  timet 

T'(Sl  x  ■  x  Sl)/Sl  a  r’(s‘  X  •••  X  s'). 

N  timet  (N—l)  timet 


Remark  3.4.2:  The  reduced  bracket  on  T*(S1  X  •••  X  51)/51  can  also  be  ob¬ 
tained  from  the  general  formula  for  the  bracket  on  ( P  X  T*G)/G  =  P  x  g*  found  in 
Krishnaprasad  and  Marsden  [34]. 

The  reduced  Hamiltonian  on  P  is  just  (3.10)  regarded  as  a  function  of  Mi,  ...,/* at 
and  Ok,j(k)>  k  =  2, We  therefore  know  that  the  Euler-Lagrange  equations 
(2.3)  are  equivalent  to  F  =  {P,  H}  for  the  reduced  bracket  (4.12). 

Reduced  Hamiltonian 

The  reduced  Hamiltonian  is  given  by 

H  =  5£Tj~‘(£)m  =  (4-17> 
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where  the  conjugate  momentum  vector  /£  =  [mi, • •  • »  /xjy]  and  6_  =  [£*,./(*),  k  =  2,  •  •  • ,  JVj 


The  momentum  map  J  is  given  by 

^(m)  =  Mi  +  +  My* 

and  the  conjugate  momenta  /x,-’s  are  given  by 

M  =  J(£)w,  (4.18) 

with  w  and  J  given  by  (2.27)  and  (2.28)  respectively. 


We  now  take  a  moment  to  review  our  analysis  and  put  the  results  obtained  so 
far  into  perspective.  The  reduced  Poisson  space  P  is  coordinatized  by  0t,j(k)i  k  = 
2, . . . ,  N  and  Mtt  k  =  1, . . . ,  N]  topologically  P  =  S1  x  •  •  •  X  S1  xTR.N.  The  Poisson 

S_.  — —  —s^,  * 

(JV— 1)  timij 

structure  on  this  manifold  is  given  by 


"  IV 

=  E  ( 

4=3  L  \ 


\( df 

df  \ 

9g 

(  dg 

_dg_) 

df 

\dfij(k) 

dUk) 

dOk,j{k) 

\dfij{k) 

Mk,j(k). 

(4.19) 


We  summarize  the  results  so  far  as  a  lemma. 

Lemma  3.4.2  :  The  reduced  Hamiltonian  for  the  planar  system  of  N  rigid  bodies 
connected  in  the  form  of  a  tree  structure  by  means  of  frictionless  revolute  joints  is 
given  by  : 

H  =  (4.20) 

where, 

M  =  (4-21) 

(V  as 

Ii  X 12  •••  Aitf 

M  M  AS 

Ai2  h  • • •  \*N 

•  •  •  • 

•  •  •  • 

•  •  •  • 

As  /v  As 

A  in  XiN  •••  /y 
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[n  „  „  1  r  n 

E  mk($kj-Su)  cos(^,i)  +  |Em*(fe  x  6ki\.Z) 


N 


h  =  h  +  e^ii^ii2- 

;=i 

The  Legendre  transform  for  the  reduced  system  maps 


sin(0jV),(4.23) 

(4.24) 


L(e  1. • — *  H(6tl{t],k  =  2,...,N,iiu...,iik),  (4.25) 


where 


fi  =  Jw. 


(4.26) 


3.5  Hamiltonian  Dynamics 
3.5.1  Dynamics 


The  dynamics  of  a  conservative  system  of  an  arbitrary  number  of  rigid  bodies  con¬ 
nected  in  the  form  of  a  tree  structure,  is  formulated  by  mapping  the  Euler-Lagrange 
equations  into  Hamilton’s  equations  with  the  help  of  the  Legendre  transformation. 
The  resulting  equations  are  consolidated  in  the  following  theorem. 


Theorem  3.5.1  :  The  Hamiltonian  dynamics  for  the  planar,  rigid  N  body  system 
with  the  bodies  connected  to  each  other  by  revolute  joints  is  given  by  : 


Vi  = 


dH 


dH 


Vj,  $.t.  J(j)=i 

dH  dH 


for  t  =  1, . . . ,  JV, 


for  *  =  2, . . . ,  N, 


(5.1) 

(5.2) 


dVi  dnJ{i) 
where  H  is  as  given  by  Lemma  3.4.2. 

Proof  :  The  Lagrangian  dynamics  in  the  absence  of  any  external  torques  or  forces 
is  given  by 


d  (dL\  dL 

dt  ~dq  =  °’  (5-3) 
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with 


L  =  -wTJw,  (5.4) 

It 

and  £  = 

The  dynamics  in  the  Hamiltonian  setting  can  be  derived  using  the  Lagrangian  for¬ 
mulation  with  the  appropriate  transformation 

/£  =  J.w.  (5.5) 

Simple  computation  using  (5.3)  -  (5.5)  shows  that, 


d3  _.  d  fdL\ 

a  =  +  3a  =  -  . 


Mr  =  _  53  +  53  A't>wywt, 

f=i  i=t+i 


where 


T»  —  dAj, 

*  d*.  * 


Differentiating  (4.20)  with  respect  to  and  using  (5.5)  we  get, 


_LL  =  _lrj.r_rfL.j-i 

d0«,j(«)  2“  rf0»,J(»)  ~ 


1  r  dJ 

=  — u>  — - OJ 


2  rf^i,/(») 

=  -  53  AyWfcWj. 


K  I,  fc,lm  . . . 


Consider, 


V  — 
v,- ..rr=J(,) 


ag 


“53  53  A«wku;, 

Vi  J(/)=‘  Vfe<*  »*»j€Sk-i,j 

+  53  A'wU>jtWj 

Vt<l,  #.«. 

JV  _  t-1 _ 

53  -  53  aLw, 
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But  this  is  the  same  as  (5.7),  therefore, 

dL 


V  —  -  dH 
d0<  ~  vi  .hw* de>.< 


It  is  obvious  from  comparing  (5.6),  (5.9)  and  (5.3)  that 
Ui  = 

Finally, 


r  m.  _  JS-  w  - 1  n 

An  aa  A»***>iV* 

y>  *.t.j(j)=i  *>J(0 


(5.9) 


(5.10) 


dH 


or. 


w  =  J  M  =  -j— , 

—  a/* 


=  0%  ~  9j(i) 

=  Uli  —  Wj(f) 


(5.11) 


Q.E.D. 


Corollary  3.5.1  :  The  dynamics  of  a  multibody  system  evolves  over  a  reduced 
Poisson  space  P  coordinatized  by  Ok,J(k)>  &  =  2, . . . ,  N  and  pit,  k  =  1, . . .  N.  Topo¬ 
logically  P  is  S1  X  •  •  •  X  S1  xlR,  .  The  system  is  Hamiltonian  in  the  Poisson  structure 


N—l  timet 

of  P  with  the  bracket  given  by 


if  >9}  =  £ 

k—2 


\(  df  . 

-  df 1 

dg 

/  dg 

dg  \  df 

duk) 

d6k,J{k) 

\dfij(k) 

dfik)  d0kiJW 

,(5.12) 


where  /,  g  :  TR.2N  1  — +  It. 

The  corresponding  dynamics  represented  in  terms  of  the  bracket  are 

Uk  =  {l*k,H}  k  =  l 
9k,j(k)  =  {Ok,j{k),H}  k  =  2,-'-,N. 

Proof  :  Use  Lemma  3.4.2  and  Theorem  3.5.1. 


(5.13) 

(5.14) 


Corollary  3.5.2  :  The  sum  of  all  the  conjugate  momentum  variables  k  = 
1, . . . ,  N  is  equal  to  the  angular  momentum  of  the  multibody  system. 
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Proof  :  The  angular  momentum  of  the  multibody  system  along  the  Z  axis  at  the 
center  of  mass  is  given  by  : 


N 

system  angular  momentum  =  ^  [h^k  +  (rj  X  . 

*=i 

From  (2.16)  and  (2.17)  we  have 

ri  =  Y.R(S,)SkJ, 

1*1 

n  =  "tmfSkj- 

;=i 


(5.15) 


\ 


Substituting  for  R(0i),  R(0j)  and  6kil  as  given  by  (1.2),  (2.19)  and  (2.18)  respectively, 
we  have 


M  x  K)  =  EE(#Ai)  x  (m)skj) 

)= 1 1=1 

-  Z  Z  { •  ^j]  cos(^j  -  ft) 

i=i  i=i 

+  [fex^jl]sin(^-^)}w*-  (5.16) 


Therefore  from  (5.15)  and  (5.16)  we  have 


N  N  N  N 

system  angular  momentum  =  Aw*  +  Z  m*  Z  Z  {  Pm  *  4,j]  cos(0;-  -  0t) 

k=l  k= 1  j=l  1=1 


N 

=  z 

k=l 


+  [\Sk,i  X  4j|]  sin(0,-  -  0()}  ujj 


N 


Ik  +  Z  mi  II  6i.k  "2 

i=i 

N  N  N 

+  Z  mk  Z  ui  Z  {  Pm  •  hj]  cos(^y  -  0i) 

k=i  y=i  <=i  w 

+  [Pm  *  Aji]  sin(^y  ~  *i)} 

N 

=  Z  Aw* 

k=  1 

N  N  N 

+  Zm*Zwy  Z  { p*,i  •  2U  cos(0,  -  ft) 

k=l  3=1  1=1  ijij 


+  [fo  x  sin(0j  -  fl()}  . 


(5.17) 
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But  from  Lemma  3.4.2  we  have  that 


AS 

h 

Ai.j 

** 

*  *  •  Ai,at 

w  1 

N 

Yl^k  =  [1,1,...,  1] 

4=1 

Al,2 

• 

• 

A# 

*  •  •  a2(^ 

•  • 

•  • 

•  • 

» 

U)  2 

• 

• 

• 

* 

H 

*•< 

A2  ,N 

•••  In 

N 

£ 
4=  1 


N 


Ik  +  X)  ^kjUk 


i=i.  ;'/* 


Substituting  for  Aj/s  from  (2.30)  in  the  above  equation  and  noting  that  A*,/  =  A y(t 
we  get 

X)^4  =  H  Zfc  +  X]  X]  H  {  pM  •  %kj]  cos (5,  -  0i) 

k= 1  k= 1  k=lj  =  U=l,  Ijij 

+  [|£*,i  x  j  []  sin(0,  -  0j)}  wy* 

Rearranging  the  summations  we  have 

XT  —  ^2  it  +  X) XI  X)  { p*.<  •  ^*j]  cosfa  -  0t ) 

4= 1  4=1  4=1  J~  1  1=1  l*j 

+  [fe  X  ?*j|]  sin(0y  -  0/)}  .  (5.18) 


Comparing  (5.17)  and  (5.18)  we  see  that 
N 

X  Mfc  =  system  angular  momentum  about  the  Z  axis.  Q.E.D.. 

J  4=1 


3.6  Internal  and  External  Torques 

A  multibody  system  in  space  can  experience  internal  and  external  torques.  An 
internal  torque  at  a  joint  acting  on  one  of  the  bodies  produces  an  equal  and  opposite 
(reaction)  torque  on  the  body  to  which  it  is  connected  to  (at  the  joint).  Torque  due 
to  a  motor  at  a  joint,  is  an  example  of  internal  torque. 

For  a  multibody  system  if  the  angular  momentum  is  conserved,  the  dynamics 
evolve  on  a  given  symplectic  leaf  of  the  Poisson  manifold  P.  In  other  words  the 
symplectic  leaf  is  determined  by  assigning  a  particular  value  to  the  system  angular 
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momentum.  Since  internal  torques  change  the  energy  of  the  system  but  leave  the 
system  angular  momentum  invariant,  the  corresponding  motion  will  be  confined  to  a 
single  leaf. 

External  torques  on  a  multibody  system  change  the  energy  as  well  as  the  angular 
momentum  of  the  system.  The  evolution  of  dynamics  switches  from  one  symplectic 
leaf  to  another  on  the  application  of  external  pulse  torques.  Examples  of  external 
torques  are,  torques  induced  by  gas  jets,  gravity  gradient  torque  ,  solar  radiation 
pressure,  solar  wind,  etc.  Actuators  may  be  used  to  create  internal  and  external 
torques  for  control  purposes. 

We  are  interested  in  studying  the  effect  of  these  torques  on  the  dynamics  of  the 
system.  The  following  theorem  gives  the  dynamics  when  internal  and  external  torques 
are  applied  to  the  multibody  system. 


Theorem  3.6.1  :  The  dynamics  of  the  planar  multibody  system,  under  the  appli¬ 
cation  of  external  and  internal  torques  are  given  by: 

„  QH  dH  x-~k  * 

k  =  E  sr  -  M —  +  E  r,_,  -  iu  +  iy*,  (e.i) 

v>,  ....  •((,>.  dei.‘  de‘J<i)  vi, .... 

for  *  =  1, . . . ,  N 

dH  dH 

0i'JM  =  diH  ~  for*  =  (6-2) 

where  T_i  is  the  internal  (reaction)  torque  acting  at  joint  t  —  1  (To  =  0),  and  T‘xt  is 
the  external  torque  acting  on  body  t. 

Proof  :  The  Lagrangian  dynamics  in  the  presence  of  any  internal  and  external 
torques  is  given  by 

dL 


jf  (dL\ 
dt  \dg) 


dg 


=  E, 


(6.3) 


where, 


L  = 

E  =  [■?!»•••  ,En]  > 


and  Fi's  are  the  generalized  force  acting  on  body  t. 
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We  know  that 


virtual  work  =  (generalized  force)  •  (displacement). 

Now,  the  total  infinitesimal  virtual  work  done  by  the  internal  torques  acting  on  the 
system,  dWtn<  is  given  by, 

SWint  =  ^2Tk-i($9j(k)  —  99k).  (6-4) 

Similarly,  the  total  infinitesimal  virtual  work  done  by  the  external  torques  acting  on 
the  system,  9Wtxt  is 

9W,xt  =  'ETri99i.  (6.5) 

1=1 

Thus  from  (6.4)  and  (6.5),  the  total  infinitesimal  virtual  work  done  by  both  internal 
and  external  torques  6W  is, 


6W  =  SWint  +  SWtxt 

N  N 

=  ~  99k)  +  'ZTrt99i. 

k=2  •= 1 

By  definition  the  generalized  force  F,  is  given  by, 


=  E  T:~i  ~  1  +  2T‘-  (6.6) 

Vj  ,.t. 

From  Theorem  3.5.1  and  (6.3)  we  have, 


-  .  _  d  (dL\ 

dL 

+  F, 

dt  ) 

dqi 

(6.7) 

and 

dL 

vy  *.«.  »=j(y) 

dH 

dH 

(6.8) 

d9J,i 

Substituting  (6.6),  (6.8)  in  (6.7)  we  get, 


^  dH  dH  _ 

H  =  E  - ^ -  +  E  Tj-X  -  Ti _!  +  T‘st , 
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Also  from  Theorem  3.5.1  we  have, 


Mo 


dH 

dUi 


Q.E.D.. 


8H 


N. 


Remark  3.6.1  :  In  the  last  result  some  special  things  happen  because  we  are  in  the 
plane  and  deal  with  pure  torques.  As  a  consequence  the  torques  seen  by  an  observer 
fixed  on  to  a  body  is  equal  to  the  torque  seen  by  an  observer  fixed  in  the  inertial 
frame. 
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CHAPTER 

_ FOUR 

FEW-BODY  EXAMPLES 


In  this  chapter  we  consider  some  interesting  examples  of  multibody  systems  in 
space.  We  begin  by  specializing  our  results  of  Chapter  three  to  the  case  of  two- 
body  and  three-body  systems.  Applications  of  our  formalism  to  terrestrial  multibody 
systems,  in  particular  to  ground  based  robot  manipulators,  is  discussed.  We  illustrate 
this  extension  procedure  for  a  planar  two-link  manipulator  on  ground.  Dynamics  of  a 
system  of  planar  N  bodies  in  space  connected  in  the  form  of  an  open  chain  is  given. 
The  notation  used  in  this  chapter  is  the  same  as  in  the  general  tree  connected  N  body 
problem  (Chapter  three). 

4.1  Planar  Two-Body  System 

A  system  of  two  bodies  in  space,  connected  together  by  a  one  degree  of  freedom 
joint  occurs  in  many  contexts.  Consider  a  two-body  system  where  the  bodies  are  of 
comparable  mass  and  inertia,  relative  to  each  other.  One  of  the  bodies  may  represent 
a  large  space  based  sensor.  Rapid  reorientation  of  the  sensor  from  one  stationary 
position  to  another  may  be  desirable.  Such  a  “step-stare"  maneuver  would  require 
thorough  knowledge  of  the  dynamics  and  thus  a  good  mathematical  model  to  formu¬ 
late  the  necessary  control  [46]. 

We  start  with  a  planar  two-body  system  as  shown  in  Figure  4.1.  The  bodies  axe 
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Figure  4.1:  Planar  two  body  system 

connected  together  by  a  frictionless  revolute  joint.  The  system  is  conservative  with 
no  external  or  internal  torques  acting  on  it.  The  centers  of  mass  of  the  bodies  are 
chosen  to  be  the  origins  of  the  respective  local  frames  of  reference.  Without  loss  of 
generality,  the  x-axes  of  the  local  frames  of  reference  of  the  bodies  is  chosen  to  be 
parallel  to  the  line  joining  the  body  center  of  mass  and  O^. 

The  Lagrangian  dynamical  equations  are  as  given  below. 


where, 


Ui  =  +  Au(0j  ,i)vl) 

Wj  =  — 7(Aw*  +  Ai2(^2,l)w|) 

9i  =  Ui 

0%  =  U>2 


(A/a  —  csAia(02,i)) 


(1.1) 
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h  =  ( h  +  ctff), 

/j  =  (72  +  €<i|)  * 


Ai2(*2.1) 

Al2(^2,l) 


m  itwq  \ 

mi+mj  J 


edid2  cos(02>i) 


(1.2) 


Remark  4.1.1  :  Here  e  is  known  as  the  ‘reduced  mass’  and  is  different  from  our 
earlier  notation  (in  Chapter  two)  that  e<  =  — for  t  =  1,2,  which  is  a  dimensionless 
quantity. 

4.1.1  Hamiltonian  Formulation 

We  now  use  Hamilton’s  principle  and  present  the  dynamics  of  the  two-body  system 
in  this  new  setting.  We  use  Lemma  3.4.2  and  Theorem  3.5.1  with  N  —  2,  02  —  [di,0] 
and  02  =  [<^2 , 0]  to  compute  the  dynamics.  The  Hamiltonian  for  this  system  is  simply 
the  total  kinetic  energy  of  the  system. 

Lemma  4.1.1  :  For  the  planar  two-body  problem 

h  =  '*•£ 

where, 

/£=  [A» 

_  A  ^12^2,1 

^12(^2tl)  72 

fix  and  /*2  are  the  conjugate  momentum  variables  and  are  related  to  the  angular 
velocities  u>i  and  w2  as  below: 

f£  =  J.[o/1,o;2]r»  (1.6) 

7x  and  I2  are  pseudo  inertia  constants  given  by  (1.2)  (  similar  to  the  ‘aguemented’ 
body  inertias  used  by  Wittenburg  [59]  ). 

Proof  :  Use  Lemma  3.4.2  with  N  =  2. 


the  Hamiltonian  is  given  by 

(1.3) 

(1.4) 

1  .  (1-5) 


45 


Theorem  4.1.1:  The  dynamics  of  a  planar  two-body  system  in  the  (02,x, Pi,Pi) 
space  is  given  by 


£i  = 

dH 

(1.7) 

902,1 

fi2  = 

dH 

802,1 

(1.8) 

bz.i  = 

8H  8H 

(1.9) 

dH2  dpi 

Proof  :  Use  Theorem  3.5.1  with  N  =  2. 

In  Chapter  5  we  will  discuss  the  numerical  simulation  results  of  the  dynamics 
of  a  planar  two-body  system  with  suitable  initial  conditions,  and  find  the  system 
equilibria.  We  analyze  the  stability  of  these  equilibria  using  the  Energy- Casimir 
method  and  show  that  this  agrees  with  the  simulation  results. 

Remark  4.1.2:  The  sum  of  conjugate  momentum  variables  (i.e.,  H\  +  ^2),  is  equal 
to  the  system  angular  momentum  -  constant  in  the  absence  of  external  torques. 

4.1.2  Poisson  Bracket 


The  phase  space  for  the  evolution  of  the  dynamics  is  on  a  Poisson  reduced  space 


P  :=  T*{Sl  XS^/S1 


(1.10) 


The  phase  space  P  is  coordinatized  by  (02,i,Pi,P2)-  The  Poisson  structure  on  the  P 
space  is  given  by: 


(i.u) 


where, 


f,9:P 


TEL. 


(1.12) 


t 
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Figure  4.2:  Dynamics  of  a  two-body  system 
Corollary  4.1.1  :  The  dynamics  rewritten  in  terms  of  the  Poisson  brackets  are: 


Mi  =  {mi,H} 

(1.13) 

M2  =  {M2,B} 

(1-14) 

rS 

**5 

*■4 

N 

'-v-' 

II 

H 

N* 

(1.15) 

Proof  :  See  Corollary  3.5.1  with  N  =  2. 

4.1.3  Control  &  Disturbance  Torques 

We  are  interested  to  study  the  effect  of  external  and  internal  torques  on  the  dy¬ 
namics  of  the  system.  These  torques  may  be  applied  as  a  form  of  control.  It  is  to  be 
noted  here  that  internal  torques  change  the  energy  of  the  system  but  the  system  an¬ 
gular  momentum  itself  is  conserved  i.e.,  the  symplectic  leaf  is  left  invariant.  External 
torques  acting  on  the  system  however  change  the  system  angular  momentum  as  well 
as  the  energy;  in  other  words  the  evolution  of  dynamics  shift  from  one  symplectic  leaf 
to  another  on  the  application  of  such  torques. 

The  dynamics  of  a  planar  two-body  system  evolves  on  the  surface  of  a  “cylinder" 
(see  Figure  4.2).  Internal  torques  leave  the  “cylinder"  invariant.  On  the  other  hand 
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the  external  torques  change  the  evolution  of  dynamics  from  one  “cylinder”  to  another. 
The  following  Corollary  illustrates  the  dynamics  when  these  torques  are  present.  It 
is  important  to  observe  the  linear  way  in  which  the  torques  appear. 

Corollary  4.1.2  :  The  dynamics  of  a  planar  two-body  system  in  the  (02,i>  Hu  fa) 
space  under  the  action  of  an  internal  control  (reaction)  torque  T  at  the  joint,  and  an 
external  torque  T(xt  on  body  1  is  given  by  : 


Ai  = 

dH  m 
—  +  T  +  T*, 

(1.16) 

dH  m 

(1.17) 

A  2  = 

<902.1  ’ 

dH  dH 

(1.18) 

02,1  = 

dn  2  dm ' 

4.2  Three-Body  Problem 

A  logical  extension  of  the  two-body  problem  is  the  three-body  problem.  A  satellite 
with  symmetrically  placed  solar  arrays  on  either  side  of  the  main  satellite  body,  may 
be  modeled  by  such  a  system.  The  results  presented  in  Section  4.1  with  respect  to  the 
Hamiltonian  formulation  can  be  extended  to  the  planar  three  body  case.  A  planar, 
three-body  system,  with  the  bodies  being  connected  in  the  form  of  a  open  chain  with 
frictionless  revolute  joints  is  considered.  It  is  assumed  that  there  are  no  external 
torques  or  forces  acting  on  the  system. 

Lemma  4.2.1  :  The  Hamiltonian  of  the  planar  three-body  problem  is  given  by 
(1.3)  with  the  conjugate  momentum  vector  f±  and  the  pseudo  Inertia  matrix  J  being 
defined  as  below: 

/£  =  (2.1) 
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9 


(2.2) 


J  = 

with 

h  = 

AT 

h  = 
= 

Ai2(^2,l)  = 

Ais(02,l  +  03,2)  = 

^23  (^3,2)  = 

A2>i(^*»j)  ~ 


«./  = 


a# 

h 

^12  (^2,l) 
Ais(^2,1  +  03, 2) 


^12(02,l)  ^is(02,l  +  #3,2) 

J2  ^23(03,2) 

^23(03,2)  X* 


II  +  (*12  +  *3l)  ||  02  ||*i 

h  +  *12  ||  S2  ||*  +  *23  ||  03  ~  OC2  IP  +  *31  ||  03  IP, 

h  +  (*23  +  *3l)  ||  S3  |p, 


C12ASa,^(02,1)  +  CS1A?„^» 

csiAs#2a(02,i  +  03,2), 

C3lAa uti-zJ02’1)  +  CsiAz»,&* 

(a  •  b)  cos (Oij)  +  |a  x  b\  sin (0,j), 


m,my 

mi+mj+ms 


*#i  *,i  =  1,2,3, 


> 


(2.3) 


a  and  6  are  any  (2  X  1)  vectors.  #2,1  and  03,2  are  the  relative  angles  between  body  2 
and  body  1,  and,  body  3  and  body  2,  respectively. 

Proof  :  Use  Lemma  3.4.2,  with  N  =  3. 


Theorem  4.2.1:  The  dynamics  of  a  planar,  three-body  system  in  the  Hamiltonian 
setting  is  given  by: 


01  = 

dH 

302,1  ’ 

(2.4) 

02  = 

a#  55" 

302,1  30s(2  ’ 

(2.5) 

03  = 

dH 

303,2  ’ 

(2.6) 

02,1  = 

a#  an 

dm  dfii  ’ 

(2.7) 

II 

N. 

d#  5if 

3//s  3/i2  ’ 

(2.8) 
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Proof  :  See  Theorem  3.5.1,  with  N  =  3. 

Remark  4.2.1:  The  sum  of  conjugate  momentum  variables  (i.e.,  Hi  +  Ht  +  fis)  is  a 
constant. 

The  phase  space  for  the  evolution  of  the  dynamics  is  a  Poisson  reduced  space  P, 
and  is  coordinated  by  ^3,2* Ml* Ms)*  The  Poisson  structure  on  the  P  space 
{., .}  for  this  case  is  given  by: 


if  i9} 

where  f,g  :P— ►  R 


(2.9) 


Corollary  4.2.1  :  The  dynamics  of  the  three-body  system  could  also  be  written  in 
terms  of  the  brackets  as  below, 


=  On,JT}  1  =  1,2,3 

02,1  =  {02, 1>#} 

03,2  =  {03,2  j 

Proof  :  See  Corollary  3.5.1,  with  N  =  3. 


(2.10) 


4.3  N-Body  Problem  (Chain) 

In  Chapter  3  we  have  given  the  Hamiltonian  formulation  of  a  planar  rigid  N-body 
system,  connected  in  the  form  of  a  tree  structure.  The  equations  for  a  open  chain  of 
rigid,  planar,  N  body  system,  can  be  realized  as  a  particular  case  of  the  tree  structure 
(see  Figure  4.3).  Note  that  J(t)  =  *  —  1  V*  =  2, . . . ,  N. 

Lemma  4.3.1  :  The  Hamiltonian  for  an  open  chain  of  planar,  N  rigid  bodies 
connected  together  by  revolute  joints  is  given  by 

H=^/£T.J~1.M,  (3.1) 
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i. 


Figure  4.3:  Planar  N-body  satellite  connected  in  the  form  of  a  open  chain 


where,  £  =  j, . . .  is  the  conjugate  momentum  vector  and  J  is  an  N  x  N 
pseudo-inertia  matrix, 


'IN 


^2  N 


In 


(3.2) 


Proof  :  See  Lemma  3.4.2. 


Theorem  4.3.1:  For  a  system  of  planar  N  rigid  bodies  connected  in  the  form  of 
a  chain,  the  dynamics  in  the  Hamiltonian  setting  under  the  application  of  internal 
torques  7J,  t  =  1, . . . N  —  1  and  external  torques  T?xt ,  *  =  1, . . .  JV,  is  given  by: 


Pi 

P‘2 


dff 

<302,1 

dH 

<303,2 


+  in  +  rr\ 

dH 

-  +  T,  -  Ti  +  t;*, 

°  02,1 
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.  dH  dH  mext 

A.  =  ^ - ~ —  +  Ti  -  iu  +  ir‘. 


PS-1  = 


dH 


dH 


+  /T1  rn  I  rr%ezt 

J-N- 1  —  -In-2  +  -ly-i, 


d$N,s-i  d$s-i  ,N- 2 

.  dH  _ 

Mat  =  - Tjv-i  +  Tn‘, 

oVn,n- i 

.  8H  dH  r  .  „ 

M  +  l,i  A  A  1 

where,  0,+i.i  is  the  angle  between  body  *  +  1  and  body  t,  for  i  =  1 ,...  ,N. 

Proof  :  Note  J(t')  =  i  —  1,  Vt  =  2,...,N  and  use  Theorem  3.6.1. 

Poisson  Structure 

The  Poisson  reduced  space  P  on  which  the  dynamics  evolve  is  coordinatized  by 
{02, u  •  •  • ,  On,n- li  Pi,  •  •  • ,  Pn)-  The  Poisson  structure  is  given  by, 


{/>?}  =  {/.s)i  -t - 1-  {/i9}w_i> 


(3.3) 


where 


if, 


(21 

df  ) 

dg 

df 

(  dg_ 

dg  \ 

\dpi 

dPi+i  J 

d0i+i,i 

d0i+l,i 

{dm 

dpi+i ) 

and  f,g:P  — ♦  R. 


(3.4) 


Corollary  4.3.1  :  The  dynamics  in  terms  of  the  Poisson  bracket  under  the  absence 
of  external  and  internal  torques,  are  given  below, 


Pi 

0»+l,« 


{pi,H}  for  *  =  1,..., N, 
{0.+i.»,  H}  for  *  =  2, ,  N, 


Proof  :  See  Corollary  3.5.1. 


(3.5) 

(3.6) 
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Figure  4.4:  N-body  (chain)  example  :  special  case 
4.3.1  N-Body  (Chain)  :  Special  Case 


Consider  a  special  case  of  the  N-body  chain  example  wherein  the  center  of  mass  of 
a  generic  body  i  is  along  the  line  joining  the  two  joints  (see  Figure  3.5).  Symbolically 
we  have, 


on  —  KiPijf\  i  —  2, . . . ,  N  2 


(3.7) 


where  «,•  6  R.  By  proper  choice  of  the  local  frames  of  reference  we  have 


(3.8) 


(3.9) 

(3.10) 


To  find  the  dynamics  of  this  system  we  first  find  the  Hamiltonian  of  the  system 
as  given  by  Lemma  4.3.1. 

Thus  we  get 


N 


/*  =  /*  +  £ 


m. 


r  2 
t>jk 


>=i 
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and 


N 

^jl  =  mk6k}-$kl  COS (0jtl), 

k=l 


-(1  -  €i)di[l,0]r  j  =  1 

'  ~  (ziij+i  Cl  +  C;  *y)  d,[l,  0]T  j  =  2, . . . ,  N  -  1  , 

-e*<fo[i,o]T  y  =  jv 

+  dy[l,0]r  1  <j  <  k-1 

Sik  +  dk[l,0}T  j  =  k 


Si:  k  +  1  <j<N 


The  dynamical  equations  are  as  given  by  Theorem  3.5.1. 


(3.11) 


(3.12) 


Symmetrical  N-Body  Case 


For  the  particular  case  of  a  planar  N-body  problem  wherein  the  bodies  have  equal 
mass,  inertia  and  other  physical  parameters,  and  with  the  center  of  mass  of  a  generic 
body  at  the  center  of  the  line  connecting  the  joints,  we  have 


m,  =  m'  Vi  =  1, . . . ,  N 
di  =  d  Vi  =  2,...,JV-l 
d\  =  dff  =  0.5  d 
k i  =  0.5  Vi  =  2, . . . ,  JV  —  1, 


the  Sicj  ’s  are  given  by 


-  (^)  *MIT  /  =  i 
.  -(>ita±2) <f[i,o]T  y  =  2,...,jv-i  , 


Si,  +  |[l,0]r  1  <  j  <  k  -  1 

Sik  +  d[l,0]r  J  =  k 


fc 
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4.4  Terrestrial  Multibody  Systems  :  An  Example 


The  dynamics  of  a  planar  multibody  system  on  ground  connected  in  the  form  of 
tree  could  be  derived  from  the  general  multibody  system  dynamics  in  the  following 
way.  Consider  a  multibody  satellite  system  where  the  main  satellite  body  is  of  infi¬ 
nite  mass.  This  could  be  easily  modeled  in  the  Hamiltonian  setting.  The  resulting 
equations  for  such  a  case  are  equivalent  to  the  dynamics  of  a  planar  ground  based 
multibody  system,  with  the  earth  playing  the  role  of  a  body  with  infinite  mass.  We 
give  the  following  procedure  for  such  a  derivation. 

The  dynamics  of  a  N  +  1  body  system  is  computed  first.  The  limiting  set  of 
dynamical  equations  as  the  body  1  mass  (and  so  its  inertia)  tends  to  infinity  gives 
the  dynamics  of  the  N  body  system  on  ground  -  with  the  absence  of  the  gravity 
term.  The  effect  of  gravity  on  the  system  can  be  modeled  separately  as  an  external 
torque.  The  procedure  is  illustrated  with  a  simple  example  of  a  planar  two-link  robot 
manipulator. 

We  derive  the  dynamical  equations  of  a  planar  two-link  manipulator  on  ground  as 
a  limiting  case  (as  described  above)  of  a  planar  three-body  satellite  system  and  verify 
it  using  an  example  in  Horn [21].  The  dynamical  equations  of  a  three  body  system  in 
space  is  given  by  (2.2)-(2.8).  As  mi  — ►  oo  so  I\  — ►  oo  we  have  from  (2.3), 


^1,2,00 

= 

m1m2 

lim  =  m2, 

mi— 00  m\  +  m2  +  m3 

Cs.l.oo 

= 

mim3 

lim  - -  =  ms, 

oo  mx  -j-  m2  -)-  ms 

ca.s.oo 

= 

m2m8  _ 

hm  -  =  0, 

"»i— 00  mi  +  m2  +  m8 

** 

Iloo 

= 

[h  +  m8  jj  a2  ||*  +  ms  ||  03  ||2] 

■fsoo 

= 

[/s  +  ms  jj  3s  || 2] , 

A* 

A12oo 

= 

m2ASai£a(*2.t)  + 

** 

^lSoo 

= 

m3Aa»,ai(03,l)> 

CO 

8 

= 

m3ASs,?8(03,2). 
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Substituting  so  in  (2.2)  we  get, 


J-1  =  Jim  J-1 

7i— *oo 


=  lim 

7i— >oo 


1 

Aoo 


h  Au  Au 
A 12  h  A23 
A23  A23  h 

0  0  0 
0  I  too  AjSoo 

0  Ajsoq  7joo 


1-1 


(4.1) 


IV  IV  Al- 

where,  =  [hoc  hoc  ~  As3oo). 

From  Lemma  4.2.1  and  (4.1),  the  Hamiltonian  for  such  a  system  can  be  calculated 


as. 


Hoc  =  Jim  7;{£TjH 

7j-*oo  * 


[7*2,7*3] 


2A 


OO 


hoo  A2300 

A2S00  hoo 


»  ■ 

7*2 

7*3 

(4.2) 


The  dynamical  equations  of  the  ground  based  planar  two-link  manipulator  in  the 
Hamiltonian  setting  computed  using  Theorem  4.2.1  are, 


7*1 

= 

dHoc 

502, 1 

7*2 

= 

dHoo 
<9*3,2  ’ 

7*3 

dHx 

503,2 

*2.1 

= 

dHco 
5t*2  ’ 

*3.2 

= 

5ifoo 

5t*s 

=  0, 


dHt 


OO 


(4.3) 

(4.4) 

(4.5) 

(4.6) 

(4.7) 


where,  H «>  is  given  by  (4.2). 


Remark  4.4.1  :  It  is  to  be  noted  here  that  =  0  could  be  interpreted  as  that  any 
movement  of  other  bodies  has  no  effect  on  the  main  satellite  body  (i.e.,  the  ground). 
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Figure  4.5:  Planar  two-link  manipulator  in  Horn’s  example 


Horn’s  Example  :  Consider  a  particular  case  of  a  planar  two-link  manipulator  as  in 
Figure  4.4.  The  links  are  in  the  form  of  two  thin,  similar  rods  of  length  /'  and  mass 
m!  with  the  centers  of  mass  of  each  link  at  the  center  of  the  rod.  The  Hamiltonian 
for  such  a  case  from  (4.2)  is  given  by, 


Hoc  = 


[a*2,  Ms] 


3 

mj  COS  (03,2) 

-1 

E  2 

mi  cos(03)2) 

m'ln 

S 

Vs 

2  [w2,Wj] 

4 

3 

|COS(03,2) 

U>2 

2 

JCOS(0S,2) 

» 

1 

3 

U>3 

p  m 

(4.8) 


We  know  that  in  the  absence  of  internal  and  external  torques  for  the  above  system, 
the  Lagrangian  L  is  equal  to  the  Hamiltonian  Hoc  (as  given  by  (4.8)).  The  Lagrangian 
dynamics  are  given  by, 


d_dL  dL 

dt  dg  dg  ’ 

where  q  =  [0i,02,03]. 

Simplifying  the  above  equation  we  get, 


#2,1 

1 

9  cos(03  2)  sin (03i2)0j  2  -  6  sin (0Si2) (02,i  +  Vs,*)2 

^3,2 

A 

{[24  +  9cos(03i2)]  Q\tl  +  [9cos(03>2)  +  6]  (02>2  +  03>2):  j  sin(03i2) 
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A  =  (9  cos’ (£3,2)  -  16), 


which  can  be  easily  verified  to  be  the  time  evolution  equations  of  a  two-link  manipu¬ 
lator  in  the  example  in  Horn  [21]. 
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CHAPTER 

_ FIVE 

RELATIVE  EQUILIBRIA  &  THEIR  STABILITY 


Stability  of  relative  equilibria  for  the  multibody  problem  is  of  interest.  Study 
of  equilibria  and  their  stability  is  essential  for  a  better  understanding  of  the  phase 
portrait  of  the  system.  Stability  of  equilibria  of  rigid  spacecraft  with  multiple  rotors 
has  been  studied  by  Krishnaprasad  and  Berenstein  [33]  in  the  Hamiltonian  setting. 
Holm,  Marsden,  Ratiu  and  Weinstein,  [18]  have  studied  the  stability  of  a  heavy  top 
under  this  setting,  using  the  energy-Casimir  method.  We  use  Arnold’s  [l]  energy- 
Casimir  method,  as  is  summarized  in  Holm,  Marsden,  Ratiu  and  Weinstein  [18]  and 
Krishnaprasad  and  Marsden  [34]  to  determine  the  equilibria  and  to  assess  their  sta¬ 
bility.  An  equivalent  alternative  to  this  method  is  to  look  for  the  critical  points  of  the 
Hamiltonian  H  restricted  to  the  symplectic  leaves  and  test  for  definiteness  of  d2H  at 
these  equilibria.  The  equilibria  axe  not  trivial  to  find  since  they  involve  the  solution 
of  a  set  of  nonlinear  algebraic  equations  with  trignometric  polynomials.  The  solution 
of  these  equations  is  a  difficult  task  and  as  such  no  general  analytical  methods  exist. 
Numerical  techniques  like  the  continuation  methods  which  in  the  past  have  been  used 
to  solve  inverse  kinematics  of  robot  manipulators  [53],  could  be  used  to  solve  these 
equations. 

Definition  :  Consider  a  dynamical  system  E: 

*  =  /U) 
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with  x  6  P  and  P  is  the  phase  space.  The  equilibria  x*  are  defined  as  the  points 
x*  E  P  such  that 

i,  =  /(x«)  =  0. 


The  energy-Casimir  technique  is  introduced  along  with  a  simple  free  rigid  body 
example.  General  formulae  to  determine  the  Casimir  functions  and  to  find  the  Hessian 
(d2H)  are  given,  for  a  general  multibody  system.  The  equilibria  for  a  two-body  system 
are  computed.  Their  stability  is  analyzed  and  the  results  confirmed  by  numerical 
simulation.  Even  though  an  upper  bound  on  the  number  of  equilibria  for  a  three- 
body  problem  can  be  found,  the  equilibria  themselves  are  difficult  to  find.  So,  a 
special  kinematic  case  of  the  three- body  problem  is  considered  and  the  corresponding 
equilibria  are  given  explicitly.  The  energy-Casimir  technique  is  applied  successfully 
to  the  case  of  N-symmetrical  bodies  wherin  the  bodies  are  rotating  with  constant 
angular  velocity  in  an  extended  position. 

5.1  Energy-Casimir  Method 


We  present  the  energy-Casimir  method  as  given  in  Holm,  Marsden,  Ratiu  and 
Weinstein  [18]  and  illustrate  it  for  a  simple  example  of  a  free  rigid  body.  The  procedure 
to  be  followed  regarding  the  method  is  outlined  in  three  steps. 


Step  1:  Consider  the  system  E  with  the  equilibrium  point  x,  whose  Liapunov  sta¬ 
bility  we  wish  to  ascertain.  Find  the  conserved  energy  H:  P—*R ,  i.e., 


dH 

dt 


=  0. 


(1.1) 


Consider  a  family  of  conserved  quantities  F  :  P—*R.  These  conserved  quantities 

% 

are  typically  Casimirs  (generated  by  symmetry  groups).  The  bracket  of  a  Casimir 
function  with  any  other  function  g ,  where  g  :  P  — ►  R,  should  be  identically  equal  to 
zero. 
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with  x  6  P  and  P  is  the  phase  space.  The  equilibria  are  defined  as  the  points 
if  €  P  such  that 

=  f(x*)  =  °- 


The  energy- Casimir  technique  is  introduced  along  with  a  simple  free  rigid  body 
example.  General  formulae  to  determine  the  Casimir  functions  and  to  find  the  Hessian 
(d2H)  are  given,  for  a  general  multibody  system.  The  equilibria  for  a  two-body  system 
are  computed.  Their  stability  is  analyzed  and  the  results  confirmed  by  numerical 
simulation.  Even  though  an  upper  bound  on  the  number  of  equilibria  for  a  three- 
body  problem  can  be  found,  the  equilibria  themselves  are  difficult  to  find.  So,  a 
special  kinematic  case  of  the  three-body  problem  is  considered  and  the  corresponding 
equilibria  are  given  explicitly.  The  energy-Casimir  technique  is  applied  successfully 
to  the  case  of  N-symmetrical  bodies  wherin  the  bodies  are  rotating  with  constant 
angular  velocity  in  an  extended  position. 


5.1  Energy-Casimir  Method 

t 


We  present  the  energy-Casimir  method  as  given  in  Holm,  Marsden,  Ratiu  and 
Weinstein  [18]  and  illustrate  it  for  a  simple  example  of  a  free  rigid  body.  The  procedure 
to  be  followed  regarding  the  method  is  outlined  in  three  steps. 

Step  1:  Consider  the  system  £  with  the  equilibrium  point  x^  whose  Liapunov  sta¬ 
bility  we  wish  to  ascertain.  Find  the  conserved  energy  H:  P  R ,  i.e., 


dH 

dt 


=  0. 


(i.i) 


Consider  a  family  of  conserved  quantities  F  :  P— *  R.  These  conserved  quantities 
are  typically  Casimirs  (generated  by  symmetry  groups).  The  bracket  of  a  Casimir 
function  with  any  other  function  g,  where  g  :  P  —>  R,  should  be  identically  equal  to 
zero. 
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Step  2:  (  First  Variation  )  Find  all  F  in  step  1  s.t.  Hr  =  H  +  F  has  a  critical 
point  at  xt  : 

d{H +  F){xJ')=Q.  (1.2) 


Step  3:  (  Second  Variation  )  The  second  derivative  (PHrixe)  *s  computed.  If 
it  is  definite  then  the  system  is  Liapunov  stable.  The  test  is  inconclusive  if  the  above 
computed  value  is  indefinite. 

We  give  here  an  application  of  the  above  method  for  the  case  of  a  free  rigid  body 
(see  Holm,  Marsden,  Ratiu,  and  Weinstein  [18]). 

Example  5.1.1  :  Free  Rigid  Body 
The  free  rigid  body  equations  are 

m  =  m  x  w,  (1.3) 

where  m,w  £  Rs,  w  is  the  angular  velocity  and  m  is  the  angular  momentum, 
m  is  related  to  w  by  m  =  .1  •  u  where  I  =  diag(/i,  I2,  h)  is  the  diagonalized 
moment  of  inertia  tensor.  The  system  is  Hamiltonian  in  the  Lie-Poisson  structure  on 
Rs  considered  as  a  dual  of  the  Lie  algebra  of  the  special  orthogonal  rotation  group 
50(3).  The  Lie-Poisson  bracket  is  given  by, 

{F,G}  =  -m  (v-Ftm)  x  vG(m))  ,  (1.4) 

where  F,  G  :  ]RS  i— ♦  ]R. 

It  can  be  easily  verified  that  for  the  bracket  the  system  is  Hamiltonian  in  the 
sense  of  m,-  =  ,  »'  =  1,2,3,  where  m  =  [mi,mj,mj]r.  The  Hamiltonian  H 

is  simply  the  kinetic  energy  of  the  system  and  is  given  by 

H  =  ^mr  •  I-1  *  m.  (1.5) 

The  Casimir  function  for  the  bracket  (1.4)  is  any  smooth  function  <f>  :  R  i— ►  R 
such  that 

c.  =  *(^).  (1.6) 
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An  easy  computation  shows  that  the  bracket  of  <j>  with  any  other  function  G  is  zero. 
Application  of  the  energy-Casimir  method  follows. 

Step  1  :  The  equilibria  of  (1.3)  occur  when  m  and  w  are  parallel  to  one  another. 
Without  loss  of  generality  the  direction  of  m  can  be  taken  as  the  x  axis.  Further  m 
could  be  normalized  as  m*  =  [1,0,0]. 

Step  2  :  (First  Variation) 


The  derivative  of 


1  1 

Hc4  -  ^ 21  1  21  +  4>{  2  ) 

(1.7) 

is 

\mr 

dHc+  •  6m  —  (u  +  mi i  (  ))  •  6m, 

it 

(1.8) 

For  the  first  variation  dHc4  to  be  zero  we  have 

*■<¥>  -  -b 

(1.9) 

Step  3  :  (Second  Variation) 

The  second  derivative  of  (1.7)  at  equilibrium  point  m*  =  [1,0,0]  from  (1.8)  and 
(1.9)  is, 

d2Hc+-6m  =  6v-6m  +  <£,(-^-)  •  |£m|2  +  {lIU  ’  6  m)2  4>"  {}-—-) 

-  ±{ITr  -  (J7r + 

=  (t,  -  z) + 

The  various  conditions  for  the  definitiveness  of  (1.10)  are  as  given  below. 

(i)  The  quadratic  form  is  positive  definite  if  and  only  if 


<t>”  >  0, 

h  >  h,  h>  h. 
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The  function 


=  -j\&\  +  »  (i-11) 

makes  d2Hc+  positive  definite  at  the  equilibrium  point  [1,0, 0].  So  any  stationary 
rotation  along  the  longest  axis  is  Liapunov  stable. 

(ii)  The  quadratic  form  is  negative  definite  if  and  only  if 

<  o, 

1\  <  Jj,  I\  <  Is. 


An  example  of  the  function  <j>  which  satisfies  the  above  condition  is 

*(M’)  =  -jrN*  -  (^-!)\  (M2) 

which  makes  d?Hc4  positive  definite  at  [1,0,0].  This  proves  that  the  rotation  around 
the  short  axis  is  (Liapunov)  stable. 

(iii)  The  quadratic  form  is  indefinite  if 

h  >  Is,  h  <  h, 

ao  we  cannot  prove  by  this  method  that  the  middle  axis  is  unstable. 

In  summary  the  motion  of  a  free  rigid  body  along  the  long  and  short  axes  is  stable. 


5.2  Stability  of  Equilibria  -  A  General  Formulation 

For  a  general  multibody  system  we  know  that  the  Hamiltonian  (which  is  nothing 
but  the  system  kinetic  energy),  as  given  by  Lemma  3.4.2  is 

B  =  (2.1) 

where  J  =  J k  =  2,...,^)  and  the  conjugate  momentum  vector  /x  = 
Also 

fx  =  Jw,  (2.2) 
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where  u>  is  the  vector  of  angular  velocities. 

Step  1  :  The  equilibria  for  the  system  can  be  determined  by  setting  the  dynamical 
equations  in  Theorem  3.5.1  to  zero,  i.e., 


^  dH  dH  „  ,  , 

—  X*  a/3  aa  ...  ~  ®  k  — 


Km  = 


V/  d9j,i  dOi,J(i) 

dH  _  dH 
dnj{k) 


(2.3) 


This  results  in 


‘ dH' 

dH 

Pm. 

« 

Pm. 

so 


and 


dH 


dQ 


*,m  I* 


point. 

From  (2.1)  and  (2.2)  we  have 


where, 


and 


—  0  k  —  2, . . . ,  N. 

(2.4) 

w0(say)  Vt',j  =  1,...JV, 

(2.5) 

=  uN  =  u0. 

(2.6) 

o 

II 

• 

£ 

(2.7) 

e  expression  is  evaluated  at  the 

equlibrium 

1 

(2.8) 

»ove  equation  we  get, 

1 

iH-rvo 

2^«wo* 

(2.9) 

U*  *  * » 1]T  » 

(2.10) 

+ - b^tf)- 

(2.11) 

Here  j. ie  is  the  system  angular  momentum. 
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The  constant  angular  velocity  with  which  the  bodies  rotate  when  the  system  is  at  any 
equilibrium  is 

2  H  . 

u0  =  - •  (2.12) 


Step  2  :  For  the  noncanonical  bracket  in  Corollary  3.5.1,  the  Casimirs  are  any 
smooth  functions  of  the  form  <f>  —  where,  fit  =  (/ii  +  •  •  •  +  Hn)-  It  is  easy 

to  check  that  the  bracket  of  <f>  with  any  other  function  g  :  t— ►  E.  is  equal  to 

zero.  Using  (2.5)  -  (2.7),  the  first  variation  of  (2.1)  at  any  equilibrium  point  can  be 
computed  as 


d(H  +  C,) 


d{H  +  cy) 

dfi 


Sfi  + 


■-[ 

=  6(±  +  53  W  (m3)  dm 

jy 

=  £  [w0  +  2 m4>{p3t)\  dm-  (2.13) 


Setting  the  first  variation  in  (2.13)  to  zero  would  result  in  a  constraint  on  the  Casimir 
function 


(2.14) 


Step  3  :  (Second  variation)  We  initially  calculate  the  Hessian  ( d7H )  and  then  find 
the  second  variation  -  d?[H  +  (7^).  Using  the  first  variation  as  in  (2.13)  we  find 


d2H  = 


daH  a*n 
a  g  a  tide 

B*H  9aH 
L  ae_9fi  atj 


Now  from  (2.1)  we  have 


d2H 

dn2 


=  J1. 


Taking  partial  derivatives  of  (2.1)  as  below  and  using  (2.6)  we  get, 


d2H 


ddk,j(k)dti 


j-1 


d9k,j{k)  “ 


r-l 


dj 


d&k,j{k) 


ru0. 


(2.15) 

(2.16) 

(2.17) 

(2.18) 
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where  the  N  X  1  vector  £  is  as  given  by  (2.10),  and, 


d2H 

1  d 

T  t-1  ^  V— 1 

11*  1  *  .I*#f 

dQk,j{k)dOi,j(i)' 

2  d$k,j[k) 

M  •»  J  M 

Tt-1 


=  — fj,T  3 

2- 


d2J 


ld0k,j{k)dOi,j(i)  dOk,J(k)  dOijii) 


J  “V. 


The  above  equation  along  with  (2.2)  and  (2.6)  gives 


d2H 

-  -'dr* 

d2J  n  dJ  j  aj 

dOk,J(k)dOi,J(i). 

2  “ 

« 

d9k,j{k)d0i,j(i)  &0k,j{k)  dOltJ{l)' 

21,(2.19) 


with  r  as  in  (2.10). 

Assembling  the  second  variation  from  the  Hessian  in  (2.15)  and,  (2.16),  (2.18)  and 
(2.19)  we  get 

J’1  +  Mr  +  -w.  fj-1  jjfi-rl 

J  0 

fc  =  2, . . . ,  N 

, .  [T-1  aJ  1r  d*H  1 

0  [  Otk, /(*)“] f  Je 

k  =  2,...,  N  Vk,l  =  2,...,N 


d2(H  +  Cj)e  = 


where 


R  = 


1  1  •••  1 

1  1  •••  1 

•  •  ,  • 

*  *  *  •  * 

•  •  •  • 

1  1  •••  1 


(2.20) 


(2.21) 


5.3  Equilibria:  Two-Body  Case 

We  refer  to  the  two-body  example  discussed  in  Chapter  4.  The  equilibria  for  the 
system  are  found  and  their  stability  analyzed.  To  search  for  the  equilibria  we  look 
directly  at  the  Hamilton’s  equations  on  P,&s  given  in  Theorem  4.1.1,  by  setting  them 
equal  to  zero. 

The  first  two  equations  turn  out  to  be 

"  dH 


Ml  =  -M2  = 


=  edidt  sln(0,,i)  =  0, 


(3.1) 
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whereas,  the  equilibrium  condition  02,i  =  0  becomes 

dH 

dpi 

SO,  I\Hl  —  e\fl2  ^  Ilt^2  ~~  tAply  or 

uj\  —  u/2  =  u/0  (constant).  (3.3) 

From  the  above  equations  (3.1  and  3.3)  we  get  the  equilibrium  conditions: 

1.  pi  =  constant, 

2.  pi  =  constant, 

3.  02,1  =  0  or  t 

We  are  interested  to  study  the  non-degenerate  equilibria,  i.e.,  the  system  angular 
velocity  u>0  ^  0.  Degenerate  equilibria  occur  when  U)0  —  0.  For  notational  simplicity 
we  will  call  the  (non-degenerate)  equilibrium  points  with  0  equal  to  0  and  ir  as  “equi¬ 
librium  points  1  and  2”  respectively.  Equilibrium  point  1  physically  corresponds  to 
the  uniform  rotation  of  the  bodies  in  an  extended  position  (Figure  5.1a).  In  a  similar 
way  the  equilibrium  point  2  corresponds  to  the  bodies  rotating  uniformly  in  a  folded 
position  (Figure  5.1b). 

5.3.1  Stability 

Energy-  Casimir  Method  Applied  to  the  Two- Body  System: 

First  Variation  :  Let  a  Casimir  C+  €F  be  chosen  as  below  : 

c,  =  *(«;), 

where,  p,  =  (/i*  +  p?).  Then  the  first  variation  of  the  Hamiltonian  H  (see  Lemma 
4.1.1)  when  computed  taking  into  account  (3.1  -  3.3)  would  result  in, 

<H  +  C*)  =  (J^  +  W')<foi+  +  [w]*9 

—  "*■  dp i  +  +  2 p,4>^j  dpj. 
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dH 

dpt 

=  0, 


(3.2) 


(3.4) 


body  1  body  2 


body  1 


(a)  stable  equilibrium 


(b)  unstable  equilibrium 


Figure  5.1:  Equilibrium  positions  for  a  Two-Body  Problem 


Applying  (1.2)  along  with  (3.4)  we  get 

d(H  +  C+)t  =  0  =*  =  - 


dH 


dni 

w0 


where  fi,  =  (^x  +  fi2), 

Second  Variation  :  The  second  variation  of  H  can  be  computed  either  from  the 
general  formula  (2.20)  gr  directly, 


d*  {H  +  C,)  = 


+  (4rfr  +  24>')  +  +  5^ 

Zfc  +  (4„:^ + 24>')  0  +  (4„;r + w)  ^ 


d7ff 

&M13 


93H 


a9H 


dni9ta<i  BnaBta,i 

Evaluating  all  the  terms  and  choosing  (2 +  4>')  =  0  we  get, 


93h_ 

aia.i1 


(3-5) 


i'(H  +  Ct),  =  j- 


I2  —  ^12(^2, 1)  0 

—  Ai2(^2,l)  ^1  0 

0  0  Au1u>2did2  cos(02>x) 


(3.6) 
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where  Ae  =  \j\I2  —  c2 A*2 . 

From  (3.6),  we  see  that  the  principal  minors  of  d2(H  +  C$)  are, 


It  t  u>ldxd2cos{0ttx) 

-7 — .  -r —  ana  - - - . 

A,’  A-  A, 

Since  I2  >  0  and  A,  >  0  we  can  see  that  the  definiteness  of  d2(H  +  C *)  depends 
on  bj\d\dt  cos(0jti).  The  following  table  gives  the  stability  of  the  non-degenerate 
equilibria  for  various  conditions. 


Item 

^2,1  —  0 

Stable 

Unstable 

d\d2  <  0 

Unstable 

Stable 

If  u0  =  0  (degenerate  equilibria)  then  H  is  not  definite  and  the  energy-Casimir 
test  fails.  In  fact  for  u„  =  0  the  number  of  equilibria  are  infinite. 


5.3.2  Simulation:  Two-Body  Case 


Numerical  simulation  of  the  dynamics  of  a  planar,  rigid,  two-body  problem,  was 
carried  out  to  supplement  our  knowledge  of  the  phase  space  and  also  to  confirm 
the  results  from  the  application  of  the  energy-Casimir  method.  The  simulation  code 
was  written  in  FORTRAN  and  run  on  a  VAX  11/785  under  UNIX1.  The  results  were 
displayed  on  a  IRIS  2400  Graphics  work-station.  The  initial  conditions  were  generated 
separately  for  a  particular  value  of  the  kinetic  energy  H  and  for  a  fixed  value  of  the 
conjugate  momentum  sum.  A  set  of  output  data  was  taken  for  sufficient  problem 
time. 

From  Lemma  4.1.1  we  have 


H 


h  ^12(^2, 1) 

Wi 

^12(^2, 1)  % 

u>2 

—  [llW>  +  2ebiCiUJiU2  cos(02,i)  +  -^2^3]  • 


(3.7) 


^NIX  is  the  trademark  of  Bell  Labs.  NJ. 
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Now  the  system  angular  momentum  fig  (equal  to  the  sum  of  the  conjugate  momentum 
variables)  is, 


v 

(Mi  +  Ms) 


1 

1 


But  from  Lemma  4.1.1  we  know  that  n  =  Ju>,  so 


1 

1 


(7X  +  edid.2  cos(0q))  <^1  +  {h  +  edid2  cos(0o))<*>2- 


i.e.. 


w2  = 


H g  —  (A  +  edid2  cos(^o))wi 


(3.8) 


(/2  +  6did2COs(0O)) 

where,  90  is  the  angle  #2,1  at  time  zero.  Substituting  for  cj2  in  (3.8)  from  (3.7),  we  get 
an  equation  of  the  form, 

aui*  +  bu>i  +  c  =  0,  (3*9) 


where, 


a  =  - 


j  {h  +  (didj  cos(fl0))2~  (ii  +  edid2  cos(g0)) 

(72  +  edid2cos(0o))2  2  (72  4*  f-d\d2  cos(0q))  1  2  ° 


b  = 


(Ii  +  edid2  cos(0o))  j  tdxd2  cos(0o ) 


c  — 


L(/2  +  edid2cos(0o)) 

h£ 


(-^2  +  td\d,2  cos(0q))J 


—  H. 


(/2  +  ed2d2  cos(0o))2 

A  set  of  physical  parameters  2  and  a  fixed  value  of  the  conjugate  momentum  sum 
H,  =  (a»i  +  A*a  )  were  assumed.  For  a  given  value  of  H,  from  (3.9)  and  (3.8),  we  get 
two  sets  of  values  for  the  [u>i,w2]t  vector  and  so  for  the  [/ii ,  /x2]r  vector.  Simulations 
were  carried  out  for  each  pair  of  these  initial  conditions  by  assuming  a  initial  value  of 
7r/2  for  02ii  i.e.,  0O  —  tt/2. 


2mi  =  125,  =  100,  d  1  =  0.8,  d2  —  0.6,  Ii  —  70,  I2  =  50,  and  fi,  =  50.  (SI  units). 
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A  coordinate  frame  consisting  of  x  =  (pi  +  n%)  cos(^j(i),  y  =  [tii  +  ta)  sin(02|i) 
and  z  =  (/ii  —  H2)  was  chosen  to  display  the  results.  The  trajectories  were  found  to 
live  on  a  cylinder  with  the  generators  parallel  to  the  z  axis.  The  radius  of  the  cylinder 
(as  expected)  was  equal  to  the  conjugate  momentum  sum.  For  a  given  value  of  the 
kinetic  energy  and  62,1  at  time  zero  the  trajectories  were  closed.  The  trajectories 
around  the  equilibrium  point  1  were  elliptical  in  nature  with  the  ellipse  being  laid 
out  on  the  surface  of  the  cylinder  (see  Figure  5.2).  The  major  axes  of  these  closed 
trajectories  (orbits)  were  along  the  circumference  of  the  cylinder.  The  opposite  ends 
of  the  major  axis  ultimately  met  at  a  saddle  point,  which  is  the  equilibrium  point 
2  (  see  Figure  5.3).  From  then  on  the  trajectories  were  closed  around  the  cylinder. 
From  Figure  5.3  we  can  see  that  there  are  two  homoclinic  orbits  from  the  unstable 
equilibrium  (equilibrium  point  2)  and  back  onto  itself. 

Application  of  internal  torque  (this  changes  the  energy  level  in  the  system)  resulted 
in  the  trajectory  jumping  from  one  orbit  to  another  but  still  remaining  on  the  surface 
of  the  ‘cylinder’  one  started  with.  Application  of  an  external  torque  would  result  in  the 
trajectory  jumping  from  one  concentric  ‘cylinder’  to  another  depending  on  whether 
angular  momentum  was  added  (bigger  radius  ‘cylinder’)  or  subtracted  (smaller  radius 
‘cylinder’). 

One  can  expect  that,  when  an  additional  third  body  is  attached  or  if  the  system 
is  forced  (by  means  of  a  joint  torque)  there  will  be  splitting  of  these  homoclinic 
orbits  resulting  in  chaotic  dynamics.  Melnikov’s  method  may  be  used  (see  Holmes 
and  Marsden[19]  )  for  such  analysis.  The  whole  process  of  these  evolving  trajectories 
suggests  that  the  existence  of  a  ‘energy  ellipsoid’  which  gradually  emerges  out,  from 
inside  the  cylinder,  and  the  trajectories  being  formed  by  the  points  where  this  ‘energy 
ellipsoid’  intersects  the  momentum  cylinder. 

5.4  Equilibria  :  Three-Body  Case 

Consider  a  planar  satellite  composed  of  three  rigid  bodies  connected  by  revolute 
joints  in  the  form  of  an  open  chain  (for  details  on  notation  refer  to  Chapter  4,  Section 
4.3).  The  local  coordinate  system  for  body  1  is  chosen  such  that  the  x-axis  is  parallel 
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Figure  5.2:  Planar  Two-Body  Simulation  -  Equilibrium  Point  1 


Figure  5.3:  Planar  Two-Body  Simulation  -  Equilibrium  Point  2 
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body  3 


Figure  5.4:  Planar  Three-Body  System 

to  the  line  joining  Oi  and  On  .  Similarly  the  coordinate  systems  for  body  2  and  body 
3  are  chosen  to  be  parallel  to  the  line  joining  02  and  On,  and,  the  line  joining  0$ 
and  0j3  respectively. 

With  such  preliminaries  taken  care  of,  we  define  the  vectors  a2,  Ss,/?2,/?3,  in  t-heir 
respective  local  coordinate  systems  to  be  :  3 2  =  [&i,0]r,  3j  =  [di,0]r,  /?2  =  [cj,0]r, 
Pi  =  32  +  [«i,  e2]r, 

The  equilibria  for  the  three-body  system  can  be  found  by  setting  the  dynamical 
equations  (in  Theorem  4.2.1)  to  be  zero.  This  results  in  the  following  equations  : 


8H 

8H 

(4.1) 

d02,i 

d9t,t  "  ’ 

hi 

=  u2  ~  uJi  =  0, 

(4.2) 

hi 

=  U>3  —  U/j  =  0. 

(4.3) 

From  the  above  equations  it  can  be  seen  that 

uji  =  w2  =  u>s  =  u0  (constant).  (4.4) 
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The  system  angular  momentum  n,,  and  the  Hamiltonian  H  are  given  by 

[yi  I*  +  2(Au(02ii)  +  A23  +  Asi(02ii  +  #3,2)) 

L«=i 

+  2(Ai2(02ii)  +  A23(03>2)  +  A3i(02>1  +  032)) 


H.  =  Wo 

H  =  -w2 
2° 


(4.5) 


=  2W°^’ 


or, 


w0  = 


2if 


(4.6) 


(4.7) 


Keeping  in  mind  Lemma  4.2.1  we  get, 


dH 


38 


2.1  J  e 


1  3 


fe7'J'1s)|, 


2  302,1 

1  t  ,  dJ  , 

“  -2"  J 


2,1 


1  T  3J 

=  —  -  w  —  u 

2  502,1 

w 


2(mi  +  m2  +  m3)  sin^2>l  +  *3-2)  +  Bl  sin^2-1)  +  Ci  cos(02,i)] 


0, 


or,  for  the  non-degenerate  case  (w0  ^  0), 


where, 


an 


l  sin(02(1 

+  #3,2)  +  Bi  sin(02ii)  +  Ci  cos(02.i)  =  0, 

(4.8) 

Ai 

=  mim3cidi, 

(4.9) 

Bi 

—  [rr*s(6i  +  ei)  +  m2&i]  miCi, 

(4.10) 

Ci 

=  m1m3cie2. 

(4.11) 

we  get, 

.5#  w2 

=  2(m,  +  m,  +  m,)  ['4‘sin(9j'1  +  "«>  +  =  0. 

(4.12) 
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where, 


B2  —  [m!(6i  +  e2)  +  m2ei]  msdi, 

C2  =  (mi  +  m2)m3die2. 


(4.13) 

(4.14) 


From  (4.8)  and  (4.12)  we  assemble  the  final  equilibrium  equations: 

A\  sin(02ii  +  03i2)  +  B\  sin(02(2)  +  C\  cos(02ii)  =  0 
A\  sin(02ii  +  #3,2)  +  B2  sin(^3i2)  +  C2  cos(03>2)  =  0 

We  now  use  a  polar  transformation  on  the  above  equations. 

Let, 


>  . 


(4.15) 


rj  =  so  sin(02,i)  =  £  [if  -  J 

so  sin(tf2(X  +  03<2)  =  i 


5  cos(tf2il)  =  |  [17  +  J] 

[t~\]  J  cos(^ti  +  «3|2)  =  |[f  +J] 

(4.16) 

Substituting  from  (4.16)  for  the  trignometric  terms  in  (4.15)  and  simplifying  we  get, 

+  +  +  +  =  0,  (4.17) 

—  —  A\  +  ( B2  +  *C2)^  +  {—B2  +  t’C2)^  =  0.  (4-18) 

Further  manipulation  of  the  above  equations  yield, 


(P7  +  197) 7  +  (pe  +  *96)76  H - f-  (pi  +  *9x)» 7  +  (Po  +  *9o)  =  0, 


(4.19) 


a  _  1  1  M  +  {~B\  +  iCi)rj 
N  ~(B2  +  *C2)  +  A\Tf 


(4.20) 


where  p,-’s  and  9,’s  are  real  constants,  and  *  is  the  imaginary  operator. 

V 

From  the  fundamental  theorem  of  algebra  we  see  that  the  polynomial  equation  (4.19) 
should  have  exactly  seven  roots  subject  to  the  condition  p^  +  xq7  ^  0.  From  (4.20) 
we  know  that  for  every  value  of  r)  we  have  two  values  of  £. 

Therefore  in  total  we  have  14  solutions  of  the  (77 ,  £)  pair 
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==>  14  solutions  of  the  {62,1,03,2)  pair 

The  obvious  question  is  whether  all  the  above  fourteen  solutions  of  (4.19)  are 
feasible,  since  we  have  not  taken  into  account  the  constraints  |  rj  \—  1  and  |  £  |=  1. 
Several  approaches  have  been  tried,  to  answer  this  question  satisfactorily  .  Since  the 
solutions  have  to  lie  on  a  unit  circle  an  analogy  to  the  stability  of  a  discrete  time 
control  system  whose  characteristic  equation  is  (4.19),  can  be  considered.  This  has 
not  led  us  to  any  conclusive  result  since  the  resulting  expressions  are  horrendous. 

Another  approach  is  the  following.  If  the  ratio  of  the  magnitudes  of  the  coefficients 
of  rj°  to  T]r  is  equal  to  one,  then  all  the  solutions  have  to  lie  on  the  unit  circle,  otherwise, 
we  know  that  there  are  only  an  even  number  of  solutions  which  are  less  than  fourteen. 

Yet  another  approach  would  be  that  of  restricting  ourselves  to  a  sub-class  of 
problems  with  simple  kinematics,  solving  them  and  hope  to  gain  some  insight  about 
the  number  of  feasible  solutions. 

With  this  motivation,  the  following  special  case  of  the  three  body  system  wherein 
the  joints  and  the  centers  of  mass  of  the  bodies  lie  in  a  straight  line  was  considered 
and  solved.  It  was  found  that  there  may  be  4  or  6  solutions  of  the  {62, 1,  #3,2).  pair  for 
this  special  kinematic  case  of  the  three  body  system. 

5.5  Three-Body  System  :  Special  Kinematic  Case 

We  consider  here  a  case  of  the  three-body  system  with  a  special  kinematic  struc¬ 
ture  where  the  centers  of  mass  of  the  bodies  are  aligned  with  the  joints  in  a  straight  line 
when  the  bodies  are  in  a  stretched  out  position.  It  follows  that /?3—a2  =  [ei,0]T. 

So  from  (4.11)  and  (4.14) 


t2  —  0  V  C\  —  Cj  —  0. 


Thus  (4.15)  reduces  to, 


A\  sin(02,i  +  #3,2)  +  B\  sin(02,i)  = 
Ai  sin(02,i  +  #3, 2)  +  B2  sin(0s>2)  = 
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0, 

0, 


(5.1) 

(5.2) 


with 


Ai  =  cidxmitns, 

(5.3) 

Bi  =  [(&i  +  ei)m3  +  &im2]cimi, 

(5.4) 

B2  =  ((hi  +  «i)mi  +  «im2]<iim3. 

(5.5) 

Subtracting  (5.1)  and  (5.2)  we  get 

sin(03i2)  =  /csin(02>i), 

(5.6) 

where 

B\ 

K  b2 

• 

(5.7) 

Expanding  (5.1)  and  substituting  (5.6),  we  get 

Aisin(52,i)  (cos(53i2)  +  <ccos(02il)  +  t]  =  0, 

(5.8) 

where 

B i 

r  M 

“  • 

(5.9) 

Consequently  from  (5.6)  and  (5.8)  we  have 

sin(02ii)  =  0  and 

sin(fl3l2)  =  0, 

(5.10) 

or, 

sin(03,2)  =  /csin(02ii), 

(5.11) 

cos(03|2)  +  /ccos(02iJ 

)  +  r  =  0. 

(5.12) 

It  is  obvious  from  considering  (5.10)  that  the  following  four  roots 

of  the  {81,1,63,2} 

pair  could  be  readily  identified  : 

{0,0} 

{0,7r} 

(5.13) 

► 

{tt,0} 

{*>*}  . 
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Figure  5.5:  Fundamental  equilibria 

We  label  these  equilibria  as  the  fundamental  equilibria.  A  stick  figure  representation 
(Figure  5.5)  helps  in  bringing  out  the  symmetrical  way  in  which  these  equilibria  occur. 

The  remaining  equilibria  for  this  system  are  computed  as  the  solutions  to  (5.11) 
and  (5.12).  Since  the  equilibrium  equations  are  nonlinear  and  parameter  dependent, 
one  needs  to  exercise  care  while  solving  them.  The  parameter  dependence  of  the 
equilibrium  solutions  can  be  summarized  by  two  sets  of  constraints  -  parameter-sign 
and  parameter-value  constraints  respectively.  It  was  found  that  two  extra  equilibria 
(other  than  the  fundamental  equilibria )  can  exist  at  a  time,  subject  to  the  existence 
of  suitable  values  of  k  and  r  satisfying  these  constraints.  The  maximum  number 
of  equilibria  for  a  general  three-body  system  (special  kinematic  case)  is  thus,  6.  For 
some  values  of  K  and  r  not  satisfying  these  constraints  and  for  the  cases  with  k  and/or 
t  being  zero  these  extra  equilibria  merge  with  the  fundamental  equilibria  to  give  a 
total  of  four  equilibria. 

5.5.1  Parameter-Sign  Constraints 

This  constraint  set  restricts  the  existence  of  values  of  {^2,i>^s,2}  pair  depending 
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jc>0  x  >0 


□  k:<0  t  <0 

n  *>°  *  <° 


K  <0  X  >0 


2k 


2.1 


Figure  5.6:  Parameter-Sign  Constraints 


on  the  signs  of  k  and  r. 
Using  (5.9)  in  (5.1)  we  get 


+  *3,l)  = 


(5.14) 


Taking  into  account  the  signs  of  K  and  r,  from  (5.11)  and  (5.14)  we  get  Figure 
5.6,  which  illustrates  the  feasible  regions  of  the  solution  pair  {^2,i»^3,2}  to  form  the 
parameter-sign  constraints. 


5.5.2  Parameter-Value  Constraints 

The  existence  of  solutions  of  (5.11)  and  (5.12)  is  also  dependent  on  the  actual 
values  of  k  and  r  (which  are  constants  for  a  given  three-body  system).  The  parameter- 
value  dependence  of  the  solutions  can  be  formulated  by  squaring  and  adding  (5.11) 
and  (5.12),  and  simplifying  to  get 

1  -  * 3  -  r3 


cos(02ii)  = 
cos(03,j)  = 


2kt 

K?  —  T*  —  1 

2r 


(5.15) 

(5.16) 
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2  2 


Figure  5.7:  Parameter- Value  Constraints 
so, 

I  _  k*  _  T  t 

1  <  2kt  <  *’  <517> 

1  <  ~  1  <  X  ■  (5.18) 

These  equations  could  be  represented  in  the  form  of  a  graph  as  in  Figure  5.7.  The 
graph  has  been  drawn  for  k'  >  0  and  r 1  >  0,  where 


5.5.3  Local  Frames  of  Reference 

It  is  necessary  to  choose  a  local  frame  of  reference  for  each  of  the  bodies  in  order 
to  parameterize  the  system  and  study  the  system  equilibria.  Refer  to  Figure  5.8. 
Proper  choice  of  the  local  frames  of  reference  for  bodies  1  and  3  results  in  the  vectors 
02  =  [ci,0]r  and  as  =  [<fi,0]T,  where  both  Ci  and  di  are  positive.  In  general,  the 
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Figure  5.8:  Reference  Configuration 

local  frame  of  reference  of  body  2  could  be  chosen  in  such  a  way  that  /?3  —  S2  = 
[ci,ej]T  =  [ei,0]r,  where  ex  is  positive.  Note  that  if  a2  =  [6i,0]T,  the  kinematic 
parameter  6X  could  be  either  negative  or  positive.  The  two  cases  of  the  signs  of  6x 
represent  whether  the  center  of  mass  of  body  2  is  inside  the  line  segment  joining  the 
joints  0\ 2  and  O23  or  outside  it.  If  any  of  the  kinematic  parameters  cx  or  d\  is  equal 
to  zero  then  the  three  body  problem  decomposes  into  a  two-body  problem  and  a  one- 
body  problem.  It  is  also  important  to  observe  that  with  this  choice  of  local  frames  of 
reference,  Ai  is  positive  (see  (5.3)). 

5.5.4  Parameter-Dependent  Equilibria 

We  now  delve  into  particular  cases  of  the  signs  of  parameters  k  and  r  and  establish 
the  solutions  to  the  equilibrium  equations.  We  constantly  refer  to  (5.3)-(5.9)  while 
formulating  the  necessary  conditions. 

In  all  the  cases  we  consider,  we  first  ascertain  that  there  exist  physically  realizable 
values  of  the  kinematic  parameters  -  cx,  b\,  t\  and  dj,  before  finding  the  actual 
solutions.  The  equilibria  are  evaluated  based  on  the  signs  of  cos(02|j),  and,  cos(03(2) 
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(see  (5.15)  and  (5.16)  ),  and  according  to  the  parameter-sign  and  parameter-value 
constraints.  The  results  are  presented  in  the  form  of  a  table  for  each  case.  The 
graphs  under  the  column  parameter-sign  constraints  have  to  be  read  with  02> 1  as  the 
X-axis  and  03(2  being  the  Y-axis  (see  Figure  5.6  for  more  details).  The  shaded  regions 
represent  the  valid  regions  of  existence  of  the  {02, i>  03,2}  pair.  In  the  column  of  the 
parameter-value  constraints,  the  regions  referred  to  are  the  regions  of  Figure  5.7. 

Given  values  of  k  and  r,  one  can  identify  the  corresponding  table  depending  on 
the  signs  of  these  parameters,  and  determine  which  region  they  belong  to  with  regard 
to  Figure  5.7.  The  two  extra  equilibria ,  if  any,  could  then  be  read  off  from  the  table. 

CASE  1  :  k  >  0  r  >  0 

For  k  and  r  to  be  greater  than  zero,  AuBi  and  B 2  should  be  greater  than  zero.By 
choice  of  the  local  frames  of  reference  we  have  from  (5.4)  and  (5.5)  : 

(6i  +  ex)m3  +  bim3  >0  so  «x  >  -  (l  +  — )  blt 

\  mj/ 

(6i  +  ei)mx  +  exm2  >0  so  et  >  -( — — — )  bu 

\mi  +  m2/ 

i.e., 

This  is  automatically  satisfied  if  b\  >  0. 

The  equilibrium  solutions  are  given  in  a  compact  form  in  Table  5.1. 

CASE  2  :  k  <  0.  r  <  0 

The  case  k  <  0  and  r  <  0  can  be  realized  if  and  only  if  B3  <  0,  and  B3>  0  (since 
A\  >  0  always). 

Simplifying  so  from  (5.4)  and  (5.5)  we  have 

-  (l  +  — )  h  >  ei  >  -  ( — — - )  &l  (5.20) 

\  m5/  \mx+m2/  v  ' 

Naturally,  the  above  equation  indicates  that  this  case  is  possible  only  if  6i  is  negative 

(since  ex  >  0). 
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TABLE  5.2  :  k  <  0  r  <  0 
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But  t\  >  0,  so  bi  has  to  be  negative.  Then  from  (5.22)  jjJj  is  greater  than  1  but  less 
than  a  fraction  -  which  is  impossible. 

So  kinematic  parameters  satisfying  K  <  0  and  r  >  0  can  never  exist. 


5.6  Equilibria  of  N-Body  (Chain)  :  Special  Case 

We  consider  the  special  case  of  the  N-body  (chain)  discussed  in  Chapter  4.3.1  with 
no  internal  or  external  torques.  The  dynamical  equations  given  by  Theorem  4.3.1  are, 


Ui 

dH 
d02,  i’ 

iii 

dH 

d0Z'2 

dH 
d$2,i 1 

in 

• 

dH 

dH 

dQi+i,i 

• 

« 

(in- i 

• 

8H 

dH 

d0N,N- 1 

Un 

dH 

90N,N- 

* 

-1 

dH 

?*  fori 

dtn+i 

dm 

where  m  ’s  are  the  conjugate  momentum  variables. 
Setting 


fa  =  0,  i  =  1, ...  ,N 
^i+i, »’  =  0,  i  —  1, . . .  |  N  1 

we  get  the  relative  equilibrium  conditions. 

Considering  non-degenerate  equilibrium  solutions  (u >„  ^  0)  we  have  on  simplifi¬ 
cation, 

u  =  [1,...  ,l]Tu;0, 
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where  u>0  is  the  constant  angular  velocity  with  which  all  the  bodies  rotate  when  the 
system  is  at  equilibrium.  Also 


dH 

90i+i,i 


N 


=  W»T  XT  A>k sin(0*j) 
;'=1  k=i+l 


0. 


(6.1) 


It  is  clear  from  6.1  that 


—  T17T,  1  —  I, 

n  =  0, 1,2, . . . 

gives  a  set  of  solutions. 

Further  it  is  obvious  that 


=  0,  *  =  1,...,JV-1, 


(6.2) 


i.e.,  the  case  wherein  the  bodies  are  in  a  stretched  out  position  and  rotating  with 
constant  angular  velocity  is  definitely  an  equilibrium  solution. 

We  now  study  the  stability  of  this  equilibrium  solution  for  the  case  of  symmetrical 
bodies  (bodies  with  same  mass,  inertia  and  physical  parameters,  and,  with  the  center 
of  mass  at  the  center  of  the  line  connecting  the  joints,  see  Chapter  4.3.1)  using  energy- 
Casimir  method. 

By  suitable  choice  of  the  second  derivative  of  the  Casimir  function  i.e., 


2  ix]+"  +  #  =  0, 


we  have  from  (2.20)  the  Hessian  at  the  equilibrium  given  by  (6.2)  as 


(6.3) 


d\H  +  Ct). 


j  1 

ViM  =  2,...,JV 


(6.4) 


We  need  to  establish  the  definiteness  of  the  above  matrix  at  the  equilibrium  point 
to  prove  stability.  Since  J  is  positive  definite  we  need  to  find  the  definiteness  of  the 
remaining  part  in  (6.4),  which  on  simplification  yields 


d2H 

1  2 

— —  _ _  / .  a 

T  a2j 

.dOk,k- id0i,i-i 

2  0 

e  ** 

d8k,k-ld0l,l-l 

(6.5) 
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where 


a2J 


=  EEU0)fcr/^*,l-2 . JV. 


and  A ’s  is  given  in  terms  of  6 ’s  (see  Chapter  4.3.1  and  Theorem  3.4.2). 

The  expressions  for  6 ’s  involve  inequalities  and  it  is  hard  to  prove  in  general  the 
definiteness  of  the  matrix  given  by  (6.5).  A  numerical  experiment  was  conducted  for 
such  systems,  with  the  number  of  bodies  varying  from  two  to  ten.  In  all  cases  the 
matrix  was  found  to  be  positive  definite.  These  results  are  given  in  Appendix  1. 

This  conclusively  proves  that  multibody  systems  consisting  of  two  and  upto  ten 
bodies  have  a  stable  equilibrium  when  all  the  bodies  are  rotating  with  a  constant 
angular  velocity,  in  an  extended  position. 
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CHAPTER 

_ SIX 

CONTROL  Si  STABILIZATION 


We  have  so  far  addressed  the  modeling  aspects  of  multibody  systems  in  space  as 
well  as  on  ground.  The  equilibria  of  these  dynamical  systems  have  been  looked  into 
and  their  stability  for  particular  cases  examined  by  means  of  rigorous  analysis  as  well 
as  simulations.  This  leads  to  the  next  logical  step  which  is  the  study  of  control  of 
multi-body  systems. 

The  issues  involved  in  the  control  of  spacecraft  can  be  divided  into  three  categories: 
sensing,  control  and  actuation.  Sensing  of  control  and  state  variables  is  essential  for 
feedback  control  and  supervisory  purposes.  Control  of  spacecraft  could  be  done  by  ei¬ 
ther  passive  (e.g.  gravity  gradient  methods)  or  active  means  (e.g.  momentum  wheels, 
gas  jets,  etc.)  [58],  Actuation  of  control  signals  could  be  achieved  by  methods  which 
conserve  mass  and  momentum  of  the  system  or  by  the  ones  which  do  not.  Techniques 
which  make  use  of  the  former  concept  are  based  on  the  use  of  momentum  and  reac¬ 
tion  wheels,  while  the  use  of  gas  jets,  magnetic  torquers  etc.  do  not  preserve  angular 
momentum.  The  control  strategy  itself  is  dictated  by  the  goal  specifications,  system 
dynamics,  external  and  internal  disturbances,  workspace  limitations,  and  lastly  the 
available  technology.  A  certain  amount  of  analysis  and  understanding  of  the  control 
theoretic  aspects  greatly  helps  in  devising  such  strategies. 

In  this  chapter  we  prove  global  controllability  as  applicable  to  such  multibody 
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systems  connected  in  the  form  of  a  tree  structure.  We  also  prove  that  these  multibody 
systems  are  linearizable  in  the  Input/Output  (henceforth  referred  to  as  I/O)  sense 
by  means  of  a  nonlinear  state  feedback.  Differential  geometric  methods  are  made 
use  of  to  achieve  this  “exact  linearization”  and  the  I/O  map  of  the  linearized  system 
is  also  given.  Finally  we  conclude  the  chapter  with  a  discussion  of  stabilizability  of 
these  systems  and  prove  a  feedback  stabilization  theorem.  In  each  section  we  will 
go  through  examples  which  amply  illustrate  the  control  problem  being  tackled  with 
reference  to  multibody  systems  -  particularly  multibody  spacecraft  and  terrestrial 
and  space-based  robots. 

6.1  Controllability 

Consider  a  robot  arm  mounted  on  a  satellite.  The  commanded  motions  of  this 
robot  arm  induce  motions  of  the  satellite  constituting  the  robot  base.  As  a  result 
the  robot  joint  angles  that  would  normally  be  commanded  (on  a  stationary  base)  to 
produce  the  desired  robot  end-effector  motion,  would  result  in  a  missed  target.  A 
simple  example  illustrates  this.  Consider  a  space  shuttle  of  mass  =  100,000  kgs. 
with  a  maximum  payload  of  =  30,000  kgs  on  the  shuttle  remote  manipulator 
arm.  We  desire  to  move  the  payload  by  say  r  =  10  m.  in  inertial  space.  Making  use 
of  standard  kinematic  relations,  with  a  point  mass  model,  the  invariance  of  center  of 
mass  position  indicates  that  the  payload  has  moved  only 

X  +  zsl.  =  7-5m-  (1-1) 

ms 

This  would  result  in  the  end-effector  missing  the  target  by  a  large  margin.  Further  if 
the  base  of  the  shuttle  arm  is  uncontrolled,  the  final  end-effector  state  is  not  function 
of  the  terminal  joint  angles  alone,  but  is  determined  by  the  entire  joint  angle  history. 
Thus  the  control  of  these  special  kinds  of  dynamical  systems  is  made  difficult  since 
the  kinematics  and  the  dynamics  are  coupled  in  these  systems. 

The  reaction  torques  and  forces  on  the  Space  Station  due  to  the  movement  of  the 
MRMS  could  be  minimized  if  the  MRMS  presents  a  low  inertia  profile  when  traversing 
the  space  station  keel.  Similarly,  for  an  open  chain  of  bodies  in  space,  it  is  of  interest 
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to  us  to  study  how  using  articulation  control  we  could  get  this  system  of  bodies  - 
which  presents  a  larger  inertia  profile,  to  form  a  closed  chain  -  which  presents  a 
relatively  smaller  inertia  profile  for  the  dynamics.  We  may  also  want  to  conduct  this 
maneuver  in  an  optimal  way  with  respect  to  time  and  fuel. 

Similar  problems  are  also  encountered  in  other  applications  involving  robotics  in 
space.  The  Telerobotic  Work  System  (TWS)  1  as  envisioned  by  NASA  [28],  is  to  be 
used  in  conjunction  with  the  Orbital  Maneuvering  Vehicle  (OMV).  The  function  of 
the  TWS  will  be  the  servicing  of  satellites;  satellite  repair;  assembly  and  construction; 
payload  handling,  and  contingency  repair  of  spacecraft.  Imagine  the  TWS  already 
docked  with  the  satellite  to  be  repaired  or  serviced,  and  has  now  formed  an  ‘open’ 
chain  along  with  the  satellite.  The  manipulator  arms  on  the  TWS  are  now  to  be 
used  to  repair  the  satellite.  Commanded  movement  of  the  arm  produces  undesirable 
reaction  forces  on  the  remaining  (connected)  bodies  including  the  satellite.  These 
forces  may  move  the  position  of  center  of  mass  of  the  satellite  and  may  rotate  the 
body  about  its  center  of  mass,  resulting  in  the  end-effector  missing  the  target. 

Thus  the  problem  of  closing  an  open  chain  of  bodies  in  space  appears  to  be  generic. 

Rapid  reorientation  of  large  space-based  sensors  from  one  stationary  orientation 
to  another  is  a  problem  of  continuing  concern  [46].  Such  a  “step-stare”  maneuver  is 
rendered  difficult  due  to  the  interaction  of  dynamics  of  the  connected  bodies.  For 
example  consider  a  space  based  telescope  which  is  connected  to  the  main  satellite 
body  by  means  of  a  revolute  joint.  Reorientation  of  the  telescope  may  be  necessary 
while  keeping  the  antennae  on  the  main  satellite  body  oriented  towards  the  earth. 

All  the  above  examples  require  that  a  suitable  control  input  be  found  to  drive  the 
system  from  any  initial  state  to  a  given  final  state  in  finite  time.  In  other  words  we 
are  interested  in  controllability  questions. 

6.1.1  Multibody  System  Controllability 

The  Hamiltonian  dynamics  of  a  planar  multibody  system  connected  in  the  form 
of  a  tree  structure  with  joint  (internal)  torques  2*,  k  =  1, . . . ,  N  —  1  and  an  external 
1  recent  telerobotic  literature  refers  to  this  as  Flight  Telerobotic  Servicer  (FTS)  [63] 
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torque  T{xt  on  body  1  is  given  by  Theorem  3.6.1  to  be 


dH 


At  -  £  jt  +  £  +  2T-, 

Vi.  ».«.  ^(i)=i  <7PX'1  Vi.  *.«.  J(i)=l 


= 


dH  dH 


Vi,  l  l.  i(i)=<  |  ( 


£  r,_,  -  T,_i, 


for  *  =  2, . . . ,  N 


f*N 


dH 


dOj v,j(n) 

dH  dH 


—  Ts- 1, 

for  *  =  2, ... ,  N. 


dm  dfijw 

The  above  equations  could  be  also  written  as, 

t  =  /(*)  +  £U 
y  =  h(x) 

where  the  state  vector  je  and  control  (input)  vector  u  are  given  by 

x  =  [mi,--’ k  —  2,  ...,JV]  , 
u  =  [r1,...,Tw_1,T1e*‘]r, 


and 


/(*)  = 


V'  OH 

2-Vi,  $.t.  /(;)= 1  otiti 

V  an  _  an 

^vi. /(;•)=<  i  aat.i<r(0 

i=2,  . . . ,  N-l 

OH 

0iN,JiN) 

0H_  _  OH  L.  —  o  xr 

L  ««  *  ~  iv 


g  = 


Q  I  w 

Qn-i,n 

and  Q  is  a.  N  X  N  —  1  matrix  whose  elements  are  either  0,  +1,  or  — 1 
[1,0, . . .  ,0]T,  and  Ojv-ixjv  is  the  N  —  1  x  N  null  matrix. 


(1.2) 

(1.3) 

1, 

(1.4) 

(1.5) 

(1.6) 

(1.7) 

(1.8) 


(1.9) 

(1.10) 
■  UL  = 
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The  output  y  is  the  vector  of  inertial  angular  velocities  w  of  the  system  and  is 
given  by  : 


V  =  h(x) 

=  w 

-  (l.ll) 

Lemma  6.1.1  :  The  matrix  Q_  is  full  rank. 

Proof  :  Let  us  first  list  some  of  the  properties  of  the  matrix  Q 

(i)  Q  is  a  N  x  N  —  1  matrix. 

(ii)  Any  column  of  Q  contains  only  two  nonzero  elements  a  (+l)  and  a  (-1), 
the  rest  are  all  zeros. 

(iii)  Consequently  adding  all  the  rows  of  Q  results  in  a  zero  row  vector.  This 
means  that  the  rank  of  Q  is  less  than  or  equal  to  ( N  —  l). 

(iv)  Addition  of  any  r  <  N  rows  of  Q  results  in  a  nonzero  row  vector  (since 
Q  represents  a  directional  vertex  matrix  of  a  connected  graph.  For  more 
details  see  Seshu  and  Reed  [48],  Lemma  5-l(b)  and  Theorem  5-1). 

(v)  The  non-zero  elements  of  the  first  row  of  Q  are  always  1. 

Let  Q.  , . . .  be  the  rows  of  Q.  Let  c*’s  be  scalar;  the  equation 

N 

ckQk  = 

4=1 

has  the  only  non-trivial  solution  (modulo  a  scale  factor), 


ci  =  c2  =s  •  •  •  =  cs  =  1,  (1.12) 

from  properties  (iii)  and  (iv). 

Thus  only  one  independent  relation  exists  among  the  rows  of  Q.  Since  Q  has  N 
rows  the  rank  of  Q  is  N  —  1.  Q.E.D. 
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Theorem  6.1.1  :  The  multibody  system  given  by  (1.6)  is  globally  controllable.,  i.e., 
one  can  find  an  input  to  take  the  system  from  any  initial  state  in  the  state  space  to 
any  given  final  state  in  arbitrarily  small  time. 


Proof  :  Differentiating  (1.11)  we  get 

y  = 


dh t 
dx~ 
dh 


=  ^[/(s)+£tt] 

dh  .  .  dh 

=  aj/te)  + 


dx: 


Again  from  (1.11)  we  have, 


dh 

TJ-  = 

j-1 

Tj-r  »J  j-l 

•  •  • 

J-1 

Q  |  w 

Qn-i,n 


=  J-1  [«  I  it]  ■ 

Recall  that  Q  Is  full  rank  from  Lemma  6.1.1,  30  [q  I  UM  is  full  rank  and 


i,e.,  from  (1.13)  we  can  choose 


u  = 


——o  =  full  rank, 
dx~ 


dh  r,  a 

-aj/tei  +  v 


dh  ' 

-l 

(1.13) 


(1.14) 


(1.15) 


where  v&U  is  the  new  input.  The  block  diagram  (Figure  6.1)  illustrates  a  realization 
by  feedback. 

So  (1.13)  becomes 

y  =  v,  (1-16) 

i.e.,  from  (1.1 1),  the  above  equation  becomes 


U)  =  v. 


(1.17) 
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Figure  6.1:  Controllable  realization 
Let 

0  =  [0i,02,*‘*,0Ar]r,  (1.18) 

where  is  the  angle  made  by  the  local  coordinate  system  of  body  i  with  the  inertial 
coordinate  system,  i.e., 

©  =  (1.19) 

From  (1.17)  and  (1.19)  we  get 

In  ©  O  ) 

+  v  I 

Qjfjf  0^(Ar  v  In  |  v  (1.20) 

y  =  u>  | 

By  choosing  a  suitable  nonlinear  control  as  in  (1.15)  we  have  reduced  the  system 
in  (1.6)  to  a  I/O  linearized  system  given  by  (1.20),  which  is  controllable. 

This  means  we  can  go  from  any  given  initial  state  [0(O),w(O)]r  to  any  given  final 
(arbitrary)  state,  [©(*/),  w(fy)]T,  in  arbitrary  time  t/  €  R  by  carefully  choosing  the 
control  input  v. 

Consequently,  we  can  go  from  any  given  initial  state  [£(0),  p(0)jr  to  any  given 
final  (arbitrary)  state,  [<?(*/),  m(*/)]  >  in  arbitrary  time  <;6R,  since 

H  =  J  v. 
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e  =  QTe, 


the  multibody  system  given  by  (1.6)  is  controllable.  Q.E.D. 


Actually  we  can  say  something  more  than  just  proving  the  existence  of  an  input 
which  takes  the  system  in  (1.20)  from  any  initial  position  at  time  zero  to  any  final 
position  at  time  T  (arbitrary).  Let  us  represent  the  system  as 


t  =  As.  +  Bn  I 
V  =  i  )  ' 


(1.21) 


where  z  =  [©,  w]T.  The  following  remark  gives  such  an  input  v *  which  could  be  used 
as  a  form  of  an  open  loop  control. 


Remark  6.1.1  :  Consider  the  system  (1.21).  One  of  the  controls  which  takes 
the  system  from  an  initial  position  z{ 0)  =  [©.(0),u>(0)]  at  time  zero,  to  a  position 
z[T)  =  [0(T),w(T)]  at  time  T  (arbitrary)  is  given  by 

v *  =  -BTe-ATrj\  (1.22) 

where  r/*  is  any  solution  of 

W(0,T)ri  =  e{0)  -  e~ATz(T),  (1.23) 

and 

W"(°,t)  =  J*  c-ATBBTe-ATrdr,  (1.24) 

for  t  6  [0,T]  is  the  controllability  grammian  [30]. 

Remark  6.1.2  :  Orientation  of  articulated  sensors  could  be  accomplished  by  first 
choosing  a  feedback  as  in  (1-15)  to  I/O  linearize  the  system  and  then  accomplishing 
the  orientation  by  means  of  the  open  loop  control  given  by  (1.22). 
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6.2  Exact  Linearization 


Exact  linearization  as  opposed  to  linearization  around  an  operating  point  involves 
the  use  of  differential  geometric  methods  to  determine  the  nonlinear  transforma¬ 
tion/feedback  which  linearizes  the  system.  Beginning  with  the  work  of  Brockett 
[4]  there  has  been  considerable  interest  in  nonlinear  feedback  theory  to  linearize  a 
nonlinear  system  by  using  the  concept  of  feedback  equivalence  among  systems  (char¬ 
acterized  via  necessary  and  sufficient  Lie-algebraic  conditions).  Systems  for  which 
feedback  equivalence  could  be  applied  and  which  could  be  linearized  can  then  be 
taken  to  be  equivalent  to  the  linear  controllable  ones.  Jackubczyk  and  Respondek 
[27]  discovered  the  necessary  and  sufficient  conditions  for  the  (local)  feedback  equiva¬ 
lence  for  a  nonlinear  system,  under  a  large  group  of  equivalence  transformation  than 
considered  by  Brockett  [4]. 

Transformation  of  a  nonlinear  system  to  an  equivalent  linear  system  has  been 
done  both  in  the  single  input  case  by  Su  [52],  and  in  the  multi-input  case  by  Hunt, 
Su,  and  Meyer  [24].  This  theory  provided  a  formal  mathematical  framework  for  the 
successful  implementation  of  an  automatic  Sight  controller  design  for  the  vertical 
and  short  takeoff  (VSTOL)  aircraft.  The  nonlinear  plant  could  be  controlled  by 
controlling  the  linear  system,  and  any  controller  design  could  be  done  on  the  “linear 
side”  of  the  system  (i.e.,  on  the  transformed  nonlinear  system).  An  inflight  testing  of 
an  automatic  flight  control  system  for  a  DHC-6  and  an  Automatic  Wing  Jet  STOL 
Research  aircraft  have  been  successful,  and  the  methodology  is  being  used  for  the  UH- 
1H  helicopter.  However  these  results  are  true  locally  only  i.e.,  ‘in  the  neighborhood 
of  the  origin’.  Hunt,  Su,  and  Meyer  [23]  have  also  formulated  a  procedure  for  the 
global  transformation  of  nonlinear  systems  to  linear  controllable  ones  for  the  single 
input  case. 

Isidori  and  Krener  [26]  consider  a  multi-input  multi-output  system  case  with  the 
inputs  entering  linearly,  and  achieved  the  Input/State  linearization  by  using  a  combi¬ 
nation  of  coordinate  changes  and  static  state  feedback.  The  resulting  transformation 
could  be  explicitly  solved  to  achieve  the  desired  linearization.  These  methods  are  in 
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general  known  as  Input/State  linearization. 

Input/Output  linearization  as  given  by  Isidori  [25]  for  a  class  of  nonlinear  systems 
with  the  inputs  entering  linearly  could  be  effected  either  by  a  change  of  coordinates 
or  by  a  proper  choice  of  nonlinear  state  feedback.  The  method  uses  Lie  brackets 
extensively  and  along  with  the  Volterra  series  expansion  achieves  the  desired  results. 
Here  we  use  the  Input/Output  linearization  procedure  as  delineated  in  Isidori[25] 
Chapter  5. 

The  rigid  body  dynamical  model  of  a  general  multi-body  system  in  space  could 
be  represented  as  , 

i  =  /(x)  +  gu, 

V  = 

where,  x  €  R2""1,  u  €  Rm,  y  €  It',  /  :  R2"'1  -4  R2*’1,  g  :  R2*’1  -4  R2""1  x 
Rm,  and  h  :  R2*"1  R*. 

Here  denotes  system  state  vector,  u  denotes  the  vector  of  control  variables  (internal 
torques  only)  and  y  is  the  vector  of  the  system  outputs,  m  and  l  are  the  number  of 
inputs  and  outputs  respectively.  It  is  to  be  noted  here  that  the  control  enters  linearly 
in  the  above  equations. 

The  technique  for  exact  linearization  basically  involves  three  steps.  Firstly  we 
determine  whether  or  not  the  system  is  Input/Output  linearizable  (here  onwards  re¬ 
ferred  to  as  I/O  linearizable  or  simply  linearizable)  by  a  possibly  nonlinear  coordinate 
transformation.  The  second  step  is  a  test  to  determine  whether  or  not  the  system  is 
I/O  linearizable  by  the  use  of  a  feedback  of  the  form, 

J i  =  <*(x)  +/3(x)n, 

where,  a  :  R2^-1  — ►  Rm  and  /?  :  R2^-1  — ♦  Rm  x  Rm  and  t;  6  Rm  is  the  new  input. 

If  this  test  is  positive,  we  go  through  the  final  step  wherein  we  determine  the 
exact  nonlinear  feedback  necessary  to  make  the  system  I/O  linearizable.  The  resulting 
nonlinear  representation  of  this  I/O  linearized  system  could  also  be  found.  It  is  to  be 
noted  here  that  the  technique  is  to  be  used  locally.  The  submanifold  over  which  this 
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method  works  successfully  could  be  found  in  the  second  step.  A  major  difference  with 
regard  to  the  model  between  this  section  and  the  previous  one  is  that  we  consider 
only  internal  torques  here. 

6.2.1  I/O  Linearization  of  the  Multibody  System  : 


The  dynamical  equations  for  the  multibody  system  with  only  internal  torques  at 
the  joints  and  no  external  torques ,  could  be  represented  as  (see  Lemma  3.4.2  and 
Theorem  3.6.1), 


l«C  IH- 

II  II 

f(x)  +  gu 
h(z)  ’ 

(2.2) 

where 

£  —  [a*i»  M2  *  *  *  ^2,i»  *  *  * » *  -  ■ ,  , 

(2.3) 

u  =  [Ti,2j  •  •  ‘Tn- 

-.f . 

(2.4) 

f[x)  is  given  by  (1.9)  and 

£  = 

2  ' 

Q-N-1,N-1 

(2.5) 

k(i)  =  [ 

k  =  2, . . .  2V] 

-  1 

Qjf-l ,N  1  Itf-l]  £  ) 

(2.6) 

where  O^j  is  a  (k  x  j)  zeromatrix,  and  Q  is  a  (N  x  N  —  1)  matrix  with  elements  0, 
-1  or  1.  Note  that  by  Lemma  6.1.1  the  matrix  Q  is  full  rank. 

LINEARIZATION : 

We  follow  a  method  outlined  by  Isidori  [25]  Chapter  5,  for  the  design  of  a  control 
system  to  re-orient  the  spacecraft  under  consideration.  We  give  here  the  following 
steps. 

STEP  I:  Check  whether  the  system  given  by  (2.2)  is  linearizable  by  a  change  of 
coordinates. 
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First  form 


T*  =  LgLjh  = 


rA  T<m 


T,i 


Im 


k  =  0, 


(2.7) 


where 


t>  =  Ltii>,K  i  =  I,-,; 

j  =  !,•••,»», 


(2.8) 


and  gi  and  hi  are  the  element  of  g  and  h(x)  respectively,  m  is  the  number  of  inputs 
(equal  to  N  —  1)  and  l  is  the  number  of  outputs  (  equal  to  N  —  1).  The  Lie  derivative 
Ljg  is  defined  as 


Ljg  = 


dg 


ld£} 


/, 


(2.9) 


with  the  following  properties  easily  proved  by  direct  computation, 

\dLrh 1 


LgLfh  = 


9 


£/»  =  k>l. 

If  T /:  =  0,  •  •  • ,  oo  is  independent  of  £  then  we  can  find  a  (possibly  nonlinear) 
transformation  to  render  the  system  I/O  linearizable. 

From  (2.7)  we  calculate  To  to  be 


To  = 


Lgh 

dh 


dz] 

=  [£y-i,jv  |  In-i] 

—  On-  i,n-i  » 

where  is  a  (N  —  1  X  N  —  l)  zeromatrix. 

Similarly  we  know  from  (2.7)  that 


Q 

Q-N-l.N-l 


(2.10) 


Tx  =  LgL/h 


(2.11) 
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But  L/h  can  be  calculated  from  (2.9  and  (2.6)  to  be, 

'dh' 


L/h  = 


dx 


fU) 


Qn-i,n  I  In-i]  •  /(*) 


an  _  an 

3^2  Bp i 


an 


=  Q 


an 

&MJ(k) 

an  _  an 
L  Ons  3hn,j(n) 

TdH 


—  dfi 


Since  from  Lemma  3.4.2 


3H 

df£ 


u 

(wi,. .  •  ,u»Ar]r  , 


we  have  from  (2.12), 

Lfh(x)  =  Qru; 

U>2  —  Wj 

—  W*  -  Wj(jfe) 

u)n  — 

=  Qt  J“V. 


Substituting  (2.13)  in  (2.11)  we  get 


Tx 


L,  (Qt J'V) 


9 


(2.12) 


(2.13) 
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=  qtj-‘  I  -  •  2 

l  at  J 

=  QTJ_1(«)Q.  (2.14) 

Since  Ti(x)  is  dependent  on  the  state  variables,  there  exists  no  transformation  to  I/O 
linearize  the  system. 

STEP  II  :  Now  we  attempt  to  I/O  linearize  the  system  using  nonlinear  feedback 
of  the  form  u(x)  =  a(x)  +  P(x)v  where  /3(x)  is  invertible  and  v  is  the  new  input. 

To  this  end  we  construct  a  formal  power  series 

T(s,x)  =  f^Tk(x)s-k~\  (2.15) 

4= 0 

where  T*  is  given  by  (2.7). 

Fact  6.2.1  :  There  exists  a  solution  to  the  I/O  linearization  problem  if  and  only  if 
T(s,x)  is  separable  i.e., 

T(a,x)  =  K(a)R(«,x),  (2.16) 

where 

K(s)  =  f^Kksk~\  (2.17) 

4=0 

and 

OO 

R(s,x)  =  R-!  +  ]TR*(s)s~*'\  (2.18) 

4=0 

and,  with  the  inverse  of  R_i  existing. 

Proof  :  See  Isidori  [25]  Chapter  5,  Theorem  (l.ll). 

NOTE:  Kv  is  a  matrix  of  real  numbers. 

Lemma  6.2.1  :  The  multibody  system  is  I/O  linearizable  by  feedback  of  the  form 
u  =  a(x)  +  /3(x)v. 

Proof  :  We  know  from  (2.7)  and  (2.15) 

T(s,x)  =  +  [Qr  *  J1  •  Q]  s  2  +  ^2  s-*-1 

4=2 

=  [lN.lS-2]  .  \qt  •  J"1  •  Q  +  £  [L„L)h]  (2.19) 

L  4=2 
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comparing  the  above  equation  with  (2.16)  we  get 


K(a)  =  Ojv-i.jy-i**  1  +  Ijv-i«  2  +  Qjr-i,N-i  •  a  3  +  •••  (2.20) 

i.e.,  Ko  —  K2  =  •  •  •  =  Koo  ==  Qln-i,n-i  and  Kj  =  —  In-i>  and 

R (s,x)  =  Qr  •  J-1  •  Q  +  £  [LgL)h\  s~k+1t  (2.21) 

k=2 

i.e., 

R-i  =  QT-J~l-Q.  (2.22) 

Thus  we  can  find  a  feedback  of  the  form 

u(x)  =  ot(x)  +  0(x)v, 

if  and  only  if  R_i  is  invertible  <=>■  R_i  is  full  rank. 

But  from  (2.21)  R_i  is  full  rank  <=>  J-1  and  Q  are  full  rank. 

It  is  clear  from  Lemma  6.1.1  that  Q  is  full  column  rank.  Since  J  is  a  pseudo  inertia 
matrix  its  inverse  always  exists,  i.e., 

==>  is  full  rank. 

So  there  exists  a  feedback  which  I/O  linearizes  the  multibody  system  around  any 
point  x 0.  The  subspace  over  which  this  is  true  is  given  by  the  subspace  over  R-i 
exists.  Since  R_i  exists  everywhere  the  system  is  I/O  linearizable  at  any  point  in  the 
phase  space. 

STEP  HI  :  To  find  &(x)  and  0{x)  using  the  structure  algorithm. 

We  apply  the  structure  algorithm  as  in  Isidori  [25]  Chapter  5,  to  find  the  expression 
for  a(x)  and  P(x).  . 

STEP  m.l  :  From  (2.10)  we  know  that 


—  Lgh  —  Oat-i, jv-i, 

(2.23) 

ViT0  = 

(2.24) 
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where 


and  Pi  is  a  matrix  with  no  rows.  K\  =  Iy_i  , 

Si  =  0  =  r0, 

=  does  not  exist, 

7l  =  K\h, 

=  [^2,1  j  •••  • 


STEP  m.2:  Now 


La7 1 

LgL/Vl 


LgLfh  —  R_i. 


Since  Ri  has  rank  /  =  (N  —  1)  the  algorithm  terminates  setting 


62  =  ri  =  (iV-1), 

P.2  = 

if*  =  does  not  exist, 

72  =  Z2-&/7l, 

=  L/A. 

STEP  in.3:  Linearizing  feedback 

With  the  introduction  of  the  feedback  input  u  =  a(x)  4-  /?(x)v 
equations  (2.2)  becomes 

i  =  /(*)  +  £  |a(x)  +  /?(k)£] 

=  /(*)  +  £a(x)  +  g/3(x)v, 


i.e. 


i  =  /(*)  +  <?Uk 

' 

y  =  M*) 


(2.25) 


(2.26) 


(2.27) 


(2.28) 


(2.29) 


(2.30) 

the  dynamical 


(2.31) 
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where 


/(*)  =  /(x)  +  £a(x),  (2.32) 

g{x)  =  gfi  (x) ,  (2.33) 

/(x),  /i(x)  and  £  axe  given  by  (1.9),  (2.5)  and  (2.6)  respectively. 

a(x):  For  the  closed  loop  system  to  be  linearizable  in  the  I/O  sense  we  have  to  have 


LjT  (x)  —  Otf-i.jy-i, 

where  /,  <7  are  given  by  (2.32)  and  (2.33)  respectively,  and, 


rfe)  = 


71 

72 

=  Lfh, 


i.e.,  (2.34)  reduces  to 

[£,r(x)]a(x)  =  -L/r(x), 

or  using  (2.13)  and  (2.35), 

[L„L/h]  a(x)  =  -£/ [QJ“V] 
But  from  (2.14)  we  have  LgL/h  =  QrJ_1Q,  therefore, 


«(*)  =  -  [qTj  *q] 


-l 


QrJ_1 1  - 


/(*)• 


(2.34) 


(2.35) 


(2.36) 


/3(x):  The  I/O  linearizing  feedback  also  needs  to  satisfy  the  following  condition, 

£;r(x)  =  [I...JO],  (2.37) 

where  q  is  the  number  of  steps  it  took  for  the  algorithm  to  terminate  (here  q  =  2). 

[i,r(x)]/3(x)  =  [i„_,  |o] 

=  (Ir.12] 

=  In-1.  *  (2.38) 
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Figure  6.2:  I/O  linearization  of  the  planar  multibody  system 
From  (2.14)  and  (2.38)  we  have 


LINEARIZED  MODEL  : 


(2.39) 


where 


t  =  /($)  +  g{x)v, 

y  =  Mss). 


(2.40) 

(2.41) 


/(x)  =  /(x) -)•  £a(x),  (2.42) 

S(&)  =  £/&(*)•  (2.43) 

/(x),  £,  fc(x),  a(x)  and  /?(x)  are  given  by  (1.9),  (2.5),  (2.6),  (2.36)  and  (2.39)  respec¬ 
tively.  Figure  6.2  gives  a  block  diagram  representation  of  the  system. 

Remark  6.2.1  •  Since  J“V  =  =  [wi,*  •  •  Jr,  we  can  rewrite  (2.36)  as 


a(i)  =  -[2tj-‘(*)2] 

=  -  [grJ-,(2)g] 


-1 


-1 


d3 


QT j_1(£)  i 


grj_1(£)  |  -  qtj-'(9) 


T-t-l, 


do 

dj 

90 


f(s) 


f{*)i 


(2.44) 


105 


Also  from  (2.39)  we  have  /3(x)  as 


/%)  =  [QtJ \  (2.45) 

From  (2.44)  and  (2.45)  it  is  clear  that  the  knowledge  of  the  inertial  angular  veloci¬ 

ties  of  the  individual  bodies  (w)  along  with  the  relative  angles  (0)  between  the  bodies 
is  enough  to  find  the  linearizing  feedback  at  any  given  instant. 

6.2.2  Output  in  terms  of  new  input: 

We  will  compute  y  in  terms  of  y  the  new  input. 

From  Isidori  chapter  5  equation  (1.3)  [25]  we  have 

y  =  W0(t)  +  f  W [t,  r)v(r)dr,  (2.46) 

Jo 

and 

Wo(t)  =  (2.47) 

n=0  n‘ 

and 

wm  =  £  £  LjLjLyhfo)  -  - ~t-  ■  <2-48) 

r>i=0  no  1  1  nl-  n°! 

Xg  -  Point  in  state  space  around  which  the  system  is  being  linearized  in  the  I/O  sense. 


Now  from  (2.47) 

Wo  (<)  =  M*>)  +  £ 

r\ —  1 

=  h^  +  Ljh^^  +  f^LyiiXg)1—,  (2.49) 

Ljhi^x^j  ~  I/f+gah(Xg}  —  L/h^Xg)  +  Lghcx^x^.  (2.50) 

From  (2.10)  and  (2.13)  we  have 

Ljh{x0)  =  grJ”V(0)  =  Qtu{  0),  (2.51) 

x*  =  [/£(O),0(O)] 

=  [mi  (0) ,  -  -  - ,  azjv  (0) ,  ^fcl  J(*)  (0) ,  fc  =  2,-.-,jV]r,  (2.52) 

w(0)  =  J_1(i(0))  •  M(0).  (2.53) 
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Also  from  (2.13),  (2.50)  and  (2.34)  we  have 

L~h(x )  =  LjLfh{x ) 

=  LfT(x) 

=  O-N-l.N-1- 

therefore, 

L~h{x 0)  =  £Ln-i%n-i  ti  >  2,  (2.54) 

and, 

Wo(0  =  0(0)  +  QTu(0)t.  (2.55) 

From  (2.48)  we  have, 

°°  (t  _  T«o 

«'(«.’■)  =  E  L?L}Lyh(f.)K  ■  /■  -7.  (2.56) 

no=0  /  7  nx!  n0! 

Using  (2.54)  and  (2.51)  we  have 

oo  no 

W[t,r)  = 

.*,=0  '  "o! 

00  _n0 

=  E  £?  .  (2-57) 

n0=0  .  n°- 

But  from  (2.38)  [LgLfh\  0(x)  =  Ijys,  substituting  so  in  (2.57)  we  have 

oo  _n0 

W(t,r)  =  E£? («-0— ; 

no-0  7  n°! 

=  Ijv_i(t  -  r).  (2.58) 

Substituting  (2.49)  and  (2.58)  in  (2.46)  we  have 

y  =  £(0)  +  ^rw(0)  +  f  fQ  T i  y{T2)dTidT2,  (2.59) 

where 

0(0)  =  [02,i  (0)»  •  •  •  »0i,/(it)(O), . . .  ,0Ar,j(Ar)(O)j  , 

w(0)  =  [w1(0),...,wt(0),...,wiv(0)]r, 
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We  sum  the  results  up  in  the  following  theorem. 

Theorem  6.2.1  :  The  multibody  system  given  by  (2.2)  can  be  Input/Output  lin¬ 
earized  by  the  application  of  static  state  feedback  u  =  oc(x)  +  /3(x)v  where  a  and 
(3  are  state  feedbacks  given  by  (2.36)  and  (2.39)  respectively,  and  v  is  the  new  input. 
The  I/O  map  with  the  feedback  is  given  by  (2.59). 


6.3  Stabilization  -  N-Body  Problem 


In  the  case  of  a  Hamiltonian  system  a  natural  approach  to  stabilization  is  pro¬ 
vided  by  Liapunov’s  second  or  direct  method  as  noted  by  several  authors  Tsinias  & 
Kaloupsidis  [54]  ,  Marino  [40],  van  der  Schaft  [56].  The  obvious  reason  being  that 
the  internal  energy  of  the  system  H  is  a  natural  candidate  for  the  Liapunov  function, 
since  ^  =  0  along  any  system  trajectory  with  the  input  equal  to  zero.  If  H  posseses 
a  strict  local  minimum  at  the  origin  this  implies  that  the  system  is  locally  stable 
[32].  However  in  the  present  context  H  does  not  possess  a  strict  local  minimum, 
but  H  +  C  does,  where  C  is  a  Casimir  function.  In  other  words  if  we  restrict  the 
system  to  the  symplectic  leaf  of  the  Poisson  manifold  then  H  has  a  strict  local 
minimum. 

The  dynamics  of  the  multibody  system  with  control  torques  at  the  joints  (in¬ 
ternal  torques)  and  an  output  y  is  given  by  (2.2)  -  (2.6). 

Consider  a  Liapunov  function  V(x)  as  shown  below, 


V{x)  =  H  +  CV 

=  (3.1) 


where  x  =  [/£r>£T]-  Now  using  (2.2)  -  (2.6)  and  keeping  in  mind  that  C+  =  (f>  {p\) 
is  the  total  system  angular  momentum),  we  have 


V{x) 


dt 


dH 

dx 


■x  + 


(3.2) 
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But 


a *  M)  _  a*  W)  vi  .  =  1  jv 

a*  Sji,  ,J  1 .  ’ 


also 


Af 


IT  A.-  =  0, 
1=1 


(3.3) 


(3.4) 


from  Theorem  3.5.1. 

Using  (3.3)  and  (3.4)  in  (3.2)  we  have 
dH  . 


V(x)  = 


dx~ 
dH  9H 
df£  80 

tdH 


(/(*)  +  gu) 
v-  dH  dH 


i= 1  U  *.LJ(j)=i  89j,i  vi  «.^?(i)=<  ’ 


h  dQ<,m 

N 

ZT.-t 

i—2 

N 


dH  dH 


dfii  dfiJ{i) 


dH  dH 


dftj(i)  dfjti 


1=2 


(3.5) 


Consider  now  a  feedback 


Ti-i  -  Ci-iOitj(i), 


on  the  system  where,  C,-  >  0 

Substituting  so  in  (3.5)  results  in 


N 


n*)  = 


»=2 


Since  V(x)  <0  and  x  =  [/£T>£T])  V'(x)  is  a  valid  Liapunov  function  for  the 

system.  Further  V  (x)  >  0.  The  system  will  converge  to  a  maximal  invariant  set 
given  by 


S  =  [x  :  V{x)  <  C,V(x)  =  O}. 
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Now 


V(x)  =  0  ==>  0  =  0 

=S>  0  =  constant  vector. 

But 

0  =  constant  =>  u>i  =  o>2  =  •  •  •  =  u>n  = 

=>•  w'  =  constant  (since  total  angular  momentum  is  conserved). 

Since  0  is  constant  and  to  =  [oji,U2,  •  •  • ,  WAr]r  is  constant  and  n  =  J_1(0)w, 

ii  is  a  constant,  i.e., 

A  =  0, 

0  =  0 

==>  5  =  {x  :  x  =  equilibrium  points  of  system  1  with  T<  =  0} 

Since  S  contains  only  discrete  set  of  points  and  no  trajectories,  the  system 
given  by  (2.2)  -  (2.6)  with  a  feedback  T, <  =  C,0,  with  C,-  >  0  is  asymptotically 

stable. 

Thus  we  give  the  following  Theorem. 

Theorem  6.3.1  :  The  planar  multibody  system  given  by  (2.2)  -  (2.6)  with  a 
feedback  =  C',_i0,ij(,)  and  C,_i  >  0  converges  to  a  maximum  invariant  set 
given  by 

S  —  {x  :  x  =  equilibrium  points  of  the  system  with  2*  =  0}  (3.6) 
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CHAPTER 

_ SEVEN 

SYMBOLIC  COMPUTATION 


We  have  used  symbolic  and  LISP  based  computation  as  tools  for  generation,  anal¬ 
ysis  and  simulation  of  dynamics  of  multibody  systems.  This  chapter  discusses  our 
efforts  in  that  direction.  In  Section  7.1  we  present  OOPSS  a  system  designed  to 
generate,  simulate  and  animate  the  dynamics  of  planar  multibody  system.  We  have 
also  used  symbolic  computation  in  a  MACSYMA  1  environment  to  generate  the  ex¬ 
plicit  dynamical  equations  of  multibody  system  using  augmented  body  method  [59] 
and  nested-body  approach  [13].  Software  has  been  written  and  validated.  Section  7.2 
discusses  the  symbolic  generation  of  dynamical  equations. 

7.1  OOPSS  :  Object  Oriented  Planar  System  Simulator 

OOPSS  -  Object  Oriented  Planar  System  Simulator,  is  a  general  purpose  soft¬ 
ware  system  architecture  designed  to  generate  the  dynamical  equations  of  a  multi¬ 
body  system  symbolically  and  to  simulate  the  equations  so  generated  numerically. 
OOPSS'animates  the  multibody  system  by  exploiting  the  high  resolution  graphics 
and  windowing  facilities  of  a  LISP  machine  and  has  been  implemented  in  Zeta-Lisp  on 
a  Symbolics  3600  series  machine.  A  nice  user  interface  is  characteristic  of  OOPSS. 
Users  can  interactively  : 

1  MACSYMA  is  a  trademark  of  Symbolics  Inc.,  Cambridge,  Mass. 
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(i)  choose  any  kinematic  or  physical  parameter  for  the  system, 

(ii)  change  any  runtime  initial  condition  -  system  energy,  system  angular 
momentum,  time  step,  maximum  problem  time,  initial  values  of  state 
and  other  variables  (angles,  conjugate  momentum  variables), 

(iii)  select  display  parameters  for  the  graphs, 

(iv)  choose  torque  laws  and  gains. 

OOPSS  is  implemented  using  the  Object  Oriented  Programming  (OOP)  tech¬ 
nique.  A  brief  introduction  to  OOP  follows. 

Objects  are  entities  that  combine  the  properties  of  procedures  and  data,  since 
they  perform  computation  and  save  local  state  [50].  Also,  objects  could  be  linked  to 
real  world  things.  A  program  could  be  built  using  a  set  of  objects.  In  OOP  we  have 
uniform  usage  of  objects  whereas  conventional  programming  uses  separate  procedures 
and  data.  Sending  messages  between  objects  causes  action  in  OOP.  Message  sending 
is  a  form  of  indirect  procedure  call  and  supports  data  abstraction. 

General  descriptions  of  objects  in  Zeta-Lisp  are  in  the  form  of  flavors.  In  Zeta-Lisp 
a  conceptual  class  of  objects  and  their  operations  are  realized  by  the  Flavor  System, 
where  part  of  its  implementation  is  simply  a  convention  in  procedure  calling  style;  part 
is  a  powerful  language  feature,  called  Flavors,  for  defining  classes  of  abstract  objects. 
Any  particular  object  is  an  instance  of  a  flavor;  for  example  a  flavor  satellite-body  may 
have  satellite-body-l  and  satellite-body- 2  as  its  instances.  Flavors  have  inheritance 
property;  thus  if  we  build  a  flavor  using  other  flavors  then  all  the  properties  of  the 
latter  are  inherited  by  the  former.  The  variables  associated  with  a  generic  object  are 
known  as  instance  variables. 

In  a  multibody  system  for  example,  a  generic  body  i  can  be  defined  as  an  object 
with  the  following  instance  variables  :  a  vector  connecting  the  joint  ( i-1 )  to  the  center 
of  mass  of  body  vector(s)  connecting  the  joint  (i-1)  and  other  joint(s)  on  body  t, 
angle  made  by  the  body  with  respect  to  the  inertial  coordinate  system  and  shape  of 
the  body. 
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Figure  7.1:  Block  diagram  representation  of  OOPSS 

Specific  operations  could  be  associated  with  the  objects  using  methods.  One  can 
create  a  method  to  define  a  specific  operation  on  any  instance  of  a  flavor  and  attribute 
special  properties  to  it.  For  instance  one  can  define  a  method  for  a  body  which  is  the 
leaf  of  a  tree  and  so  has  only  one  joint  (i.e.,  only  one  body  attached  to  it  -  contiguous 
to  and  inboard)  where  as  a  generic  body  has  two  or  more  bodies  attached  to  it. 

Functions  are  used  to  send  messages  to  instances  of  flavors  through  the  already 
defined  methods. 

OOPSS  uses  OOP  along  with  symbolic  manipulation  to  formulate  and  simulate 
the  dynamics  of  a  planar  multibody  system.  Figure  7.1  shows  the  block  diagram 
representation  of  OOPSS.  A  detailed  description  of  each  block  is  given  in  the  following 
paragraphs. 

7.1.1  Symbolic  manipulation 

The  symbolic  equation  generator  generates  the  dynamics  of  a  planar  multibody 
system  connected  in  the  form  of  a  tree  structure  in  the  Hamiltonian  setting.  The 
formulation  and  the  notation  is  the  same  as  detailed  in  Chapter  three.  This  block  is 
implemented  in  MACSYMA.  The  input  data  for  this  block  consists  of  the  following: 
the  label  of  the  body  contiguous  and  inboard  of  the  body  under  consideration  •/(*), 
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/*  Mode:  flacsyna;  Fonts:  cptfontb  —  *  -  */ 

block! 

n :  3, 

=0t23:[1.2J, 
s0 [33 : 11,2, S3. 

JJ[23:1. 

JJC3]:2. 

for  1:1  thru  n  do  1nert[l3 : 1n[13. 

beta_t[23 :natr1x( [c_13, [83), 
beta_t[33  :natr1i<(  [b_l*e_l3.  [e_23  )* 

a3pha_t  [23:natr1*<(  [b_l3,  [83  ). 
a1pha_t[33 :natr1n( [d_l3» [83) 

);■ 

Figure  7.2:  Data  file  for  the  MACSYMA  program 

the  various  kinematic  parameters  included  in  the  vectors  a ,  and  3y,  the  physical 
parameters  like  the  mass  and  inertia  of  the  bodies;  any  external  and  internal  torques. 
A  sample  data  file  is  as  shown  in  Figure  7.2  for  a  planar  three-body  system.  A  session 
of  the  automatic  equation  generation  is  given  for  the  above  data  file  in  Appendix  2. 

7.1.2  FORTRAN  simulator 

The  FORTRAN  simulator  simulates  the  dynamical  equations  generated  by  the 
symbolic  equation  generator,  and  has  been  implemented  by  FORTRAN-77  running 
on  the  LISP  machine.  The  FORTRAN  simulator  needs  numerical  values  of  all 
parameters  to  be  in  an  input  data  file.  This  input  data  file  is  generated  by  the 
DESCRIPTOR.  The  input  data  file  contains  the  numerical  values  of  all  kinematic 
and  physical  parameters,  the  system  angular  momentum  and  system  energy  values, 
problem  time,  time  step  etc.,  associated  with  the  particular  example. 

The  state  and  related  variables  (for  example,  angular  velocities)  at  any  instance 
of  time  could  be  passed  onto  the  DISPLAY  block  by  means  of  lispfunctions  to  be 
used  for  animation  and  display  purposes.  For  example  the, 

lispfunction  displaybody  displaybody  (omegal,omega2>omegaS,theta21)thetaS2,time) 
can  be  used  as 
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displaybody(u>i ,  w2,  w3,  02i,  #32,  t) 

in  the  FORTRAN  program  to  pass  the  relevant  variables  to  the  DISPLAY  block. 
The  function  ‘ displaybody'  is  implemented  in  Zeta-Lisp  in  the  DISPLAY  package. 

The  initial  condition  for  initiating  the  simulation  is  chosen  as  follows.  An  initial 
system  energy  H,  system  angular  momentum  and  relative  angle  02,i  =  0,  are 
selected.  Since  the  energy  is  quadratic  in  the  conjugate  momenta  (see  Lemma  3.4.2) 
and  the  system  angular  momentum  fi,  =  we  get  two  solutions  for  the  (p\ ,  p-f) 

pair.  User  has  the  choice  to  select  either  the  ‘First  initial  condition’  or  the  ‘Second 
initial  condition’  (see  Figure  7.3). 

7.1.4  DESCRIPTOR 

The  DESCRIPTOR  consists  of  descriptions  of  various  flavors  to  implement  the 
display  and  the  user  interface.  It  also  contains  flavors  to  define  a  generic  body  in  the 
multibody  system.  Using  methods  we  can  attribute  special  properties  to  the  instances 
of  the  flavors. 

The  user  interface  consists  of  a  window  with  many  panes  (Figure  7.3).  Three 
frames  of  references  in  the  corresponding  window  panes:  inertial  frame  of  reference 
and  two  other  selectable  frames  (from  various  joint  and/or  body  frames),  have  been 
implemented  .  The  message  pane  where  a  part  of  runtime  data  for  the  FORTRAN 
program  is  displayed.  The  simulation  pane  draws  various  graphs  of  state  and/or  other  v 
variables  as  functions  of  time.  The  menu  pane  is  self  descriptive.  Every  item  in  this 
pane  is  mouse  sensitive  (mouse  selectable).  More  on  this  in  the  following  sections  as 
we  deal  with  particular  examples  of  planar  two-body  and  planar  three-body  systems. 
A  brief  online  HELP  facility  exists  and  information  can  be  got  by  clicking  left  using 
the  mouse  when  it  is  highlighted. 

7.1.5  DISPLAY 

The  DISPLAY  block  is  the  implementation  of  various  functions  to  drive  the 
instances  of  flavors  by  sending  messages  to  them.  DISPLAY  keeps  track  of  sending 

) 
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Figure  7.3:  OOPSS  window 


proper  messages  to  the  relevant  panes  as  and  when  it  receives  data  from  the  FOR¬ 
TRAN  block.  DISPLAY  is  characterized  by  a  ‘fu:a/u-xor’  option  which  helps  in 
*  erasing  the  display  at  time  t  and  creating  a  new  display  at  time  t  +  At.  ‘ displaybody' 
and  ‘ cleanbody'  are  functions  to  display  the  system  and  clean  the  displayed  picture 
off  the  relevant  pane  respectively. 

7.1.6  Implementation 

Two- Body  Example 

Refer  to  Chapter  four  for  notation  and  description  of  the  problem.  Figure  7.4a 
shows  the  menu  pane  for  this  example.  Clicking  left  on  the  mouse  when  system 
parameters  is  highlighted  gives  Figure  7.4b.  Clicking  left  on  Torque  Law  results  in 
Figure  7.4c.  Further  clicking  on  P-D  Torque  Law  implements  a  joint  torque  law 
(internal  torque)  -  a  proportional  sinusoidal  biased  spring  plus  derivative  controller, 
i.e.,  {Kp  sin(02il  -  0bias)  +  Kd& 2,i).  Gains  Kp  and  Kd  and  the  bias  angle  could 
be  chosen  (see  Figure  7.4d)  interactively  by  the  user. 

The  stable  equilibrium  is  displayed  in  Figure  7.5.  The  trace  shown  in  the  inertial 
frame  is  the  trace  left  by  the  joint  as  it  moves  in  space.  Simple  calculation  shows 
that  this  trace  is  indeed  a  circle  when  the  system  is  in  stable  equilibrium  position. 
Figure  7.6  displays  a  trajectory  when  the  system  is  at  a  point  very  near  the  saddle 
point.  Figure  7.7  describes  a  trajectory  for  an  arbitrary  value  of  the  system  energy.  If 
a  P-D  torque  is  introduced  in  the  system  then  the  resulting  trajectory  is  as  shown  in 
Figure  7.8  -  7.9  (no  bias)  and  Figure  7.10  -  7.11  (bias).  Notice  that  with  Kp  equal  to 
zero  and  Kj,  positive  the  system  always  goes  to  the  stable  equilibrium  and  confirms 
our  ‘stabilization  theorem’  •  Theorem  6.3.1,  i.e.,  introduction  of  a  feedback  internal 
torque  proportional  to  the  rate  of  the  relative  angle  stabilizes  the  system.  One  could 
also  interpret  this  result  as  follows:  by  introducing  this  torque  law  the  energy  in  the 
system  is  dissipated  till  the  system  goes  to  a  minimum  energy  state  which  is  the  stable 
(stretched  out)  equilibrium. 
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Three-Body  Example 


A  general  three-body  example  has  also  been  implemented  on  OOPSS.  Figure  7.12 
shows  a  general  three-body  system  wherein  the  center  of  mass  of  the  middle  body  is 
not  along  the  line  joining  the  two  joints.  The  filled-in  circles  represents  the  center  of 
mass  of  each  body  (the  first  body  represented  by  a  big  circle,  second  with  a  smaller 
circle  and  the  third  with  the  smallest  circle).  Display  frames  could  be  chosen  by 
clicking  left  on  the  ‘Choose  display  frames’  using  the  mouse.  Figure  7.13  displays 
special  kinematic  case  of  a  three-body  example  discussed  in  Chapter  5.3.  Joint  torques 
of  the  proportional  sinusoidal  bias  spring  plus  derivative  type  can  be  introduced  at 
the  joints  (see  Figure  7.14  ). 

Complex  Multibody  Examples 

Plans  are  underway  to  implement  complex  multibody  system  examples.  As  the 
complexity  of  the  examples  grow  one  may  find  oneself  limited  by  the  processing  ca¬ 
pabilities  of  the  currently  available  LISP  machines.  One  could  take  advantage  of  the 
existing  parallelism  in  these  problems  by  utilizing  the  processing  power  of  parallel 
LISP  processors.  A  Connection  Machine  2  may  just  serve  the  purpose.  Thus  dy¬ 
namics  of  complex  multibody  systems  may  be  generated  automatically,  simulated  and 
animated. 

7.2  Symbolic  Dynamical  Equation  Generation 

The  formulation  of  dynamical  equations  by  hand  is  a  tedious  process  and  of¬ 
ten  prone  to  errors.  Many  researchers  have  taken  up  the  possibility  of  computer- 
aided  methods  to  generate  these  equations.  One  way  of  categorizing  these  computer- 
oriented  methods  is  by  the  way  of  numerical  and  symbolic  programs.  Numerical  pro¬ 
grams  are  characterized  by  numerical  digital  computation  whereas  symbolic  programs 
generate  equations  and  accompanying  expressions  in  symbolic  (or  alpha-numeric) 

Manufactured  by  Thinking  Machines  Inc,  Cambridge,  Mass. 
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Figure  7.5:  Two-body  problem  :  stable  equilibrium 
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Figure  7.6:  Two-body  problem  :  unstable  equilibrium 


Simulation  Data 
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Figure  7.7:  Two-body  problem  : 
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Figure  7.8:  Two-body  problem  :  K4 
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Figure  7.10:  Two-body  problem  :  Kp  =  50,  Kp  =  50,  bias=30 


126 


Figure  7.11:  Two-body  problem  :  Kp  =  50,  Kp  =  50,  bias=90 
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Figure  7.12:  General  three-body  problem 
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Figure  7.13:  Three-body  problem  :  special  kinematic  case 
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Figure  7.14:  Three-body  problem  :  special  kinematic  case  with  joint  torques 


form  on  the  basis  of  alpha-numeric  data.  Efforts  in  this  area  have  lead  to  a  vari¬ 
ety  of  simulation  tools.  We  mention  here  a  few  such  tools.  This  list  is  by  no  means 
exhaustive. 

(1)  NBOD  :  A  program  written  exclusively  in  FORTRAN  [13],  to  simulate  the 
dynamics  of  a  multi  -  rigid  body  system,  connected  in  the  form  of  a  tree  structure 
[14],  gyrostats  and  spring  loaded  point  masses  could  also  be  taken  into  account. 

(2)  TOAD  (  Tele  Operator  Arm  Design  )  :  This  program  is  written  in  PLl/FORMAC 
for  a  multi-link  robot  manipulator  with  revolute  and  prismatic  (translational)  joints 
[51]. 

(3)  OSSAM  (  OHIO  State  Symbolic  and  Algebraic  Manipulator  )  :  It  is  written  in 
LISP  and  was  designed  to  reduce  the  use  of  memory.  It  can  be  implemented  on  a 
small  computer  such  as  PDP-11  [10]. 

(4)  EDYLMA  :  (Equation  Dynamiques  Litterales  d’  un  Manipulateur  Artticule’) 
Written  for  revolute  joints  only,  in  PLl/FORMAC  [37]. 

(5)  EG  AM  :  (Equation  Generation  of  Articulated  Mechanisms  )  Written  in  PLl  / 
FORMAC  for  both  revolute  and  prismatic  joints.  An  evolution  of  EDYLMA  [37]. 

(6)  DYMIR  :  Written  in  REDUCE,  it  can  take  into  account  revolute,  as  well  as, 
prismatic  joints  [5]. 

(7)  ARM  :  Algebraic  Robot  Modeler,  generates  forward  solution  and  Lagrangian 
dynamic  robot  model;  implemented  in  ‘C’  and  ‘LISP’  [42]. 

(8)  MESA- VERDE  :  (MEchanism,  SAtellite,  VEhicle,  and  Robot  Dynamics  Equa¬ 
tions)  Simulates  nonlinear  dynamics  and  kinematics  of  articulated  multibody  systems, 
uses  augmented  body  method.  Written  in  PASCAL  (ANSI  1983  standard). 

(9)  DYNAMAN  :  (DYnamics  and  Analysis  of  MANipulators)  Program  to  auto¬ 
matically  generate  the  dynamical  equations  and  the  Jacobian  of  a  multi-link  robot 
manipulator.  Written  in  MACSYMA  it  can  also  generate  FORTRAN  code  to  simulate 
the  dynamical  equations  so  generated  [49]. 

Symbolic  programs  in  general  are  more  efficient  in  terms  of  running  time  in  com¬ 
parison  with  numerical  programs  due  to  the  following  reason.  A  symbolically  gen¬ 
erated  expression  for  a  specific  multibody  system,  is  composed  of  non- zero  variables 
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with  fixed  addresses.  Its  numerical  evaluation  takes  little  time.  In  a  numerical  pro¬ 
gram,  large  number  of  summations  and  multiplications  of  terms  which  are  zero  for  a 
specific  multibody  system  can  cause  substantial  loss  of  computation  time.  Further¬ 
more,  computation  of  addresses  for  a  large  number  of  array  elements  with  variable 
indices  can  take  significant  amount  of  computation  time.  Complicated  multibody 
system  evaluation  based  on  symbolic  methods  can  be  faster  by  a  factor  of  ten  and 
more  as  compared  to  numerical  programs  [59]. 

Symbolic  programs  have  the  ability  to  compute  “sensitivity”  of  the  model  with 
respect  to  parameters  such  as  mass,  length  or  inertia.  This  would  be  useful  for  design 
purposes.  A  symbolic  tool  allows  us  to  simulate  some  nonlinear  control  techniques 
which  would  be  difficult  to  realize  using  numerical  methods  [5],  Many  these  programs 
are  computationally  transparent  and  so  are  easier  to  implement  and  to  debug. 

Increasing  availability  of  LISP  machines  has  lead  to  the  popularization  of  LISP 
based  symbolic  tools.  These  machines  support  a  variety  of  languages  including  MAC- 
SYMA,  PROLOG,  REDUCE,  FORTRAN-77  and  C.  The  windowing ,  object  oriented 
programming  and  the  menu  driven  programming  capabilities  make  the  programming 
environment  attractive.  Excellent  debugging  facilities  facilitates  faster  turn-around 
time.  Evaluation  of  slices  of  the  program  inside  the  editor  is  possible.  Reasonable 
graphic  capabilities  and  good  resolution  are  characteristic  of  such  machines.  A  variety 
of  fonts  help  display  of  equations  in  a  compact  form  s. 

Symbolic  tools  could  be  used  to  explore  the  phase  space  characteristics  of  dynam¬ 
ical  equations  -  like  the  computation  of  equilibria  and  also  to  study  their  stability 
properties.  They  could  also  be  used  as  advanced  (symbolic)  “calculators”  to  compute, 
simplify  and  resolve  symbolic  expressions.  In  short  a  symbolic  tool  can  be  used  as  a 
discovery  tool. 

A  number  of  symbolic  computation  languages  like  MACSYMA,  REDUCE,  PRO¬ 
LOG,  MAPLE,  and  SCRATCHPAD  H  are  available.  We  chose  to  use  MACSYMA  to 
implement  our  symbolic  computation.  It  has  a  nice  interactive  facility,  it  can  be  imple- 

3Symbolics  3600  series,  T1  explorer,  Xerox  Dandelion  are  some  examples  of  LISP  machines 
available  in  the  market. 
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mented  in  any  computer  system  which  runs  LISP,  is  easily  available  and  consequently 
widespread  in  its  usage.  We  have  used  MACSYMA  to  compute  the  equilibrium  equa¬ 
tions  of  a  planar  two-body  and  planar  three-body  systems.  The  stability  of  some  of 
these  equilibria  has  been  computed  using  energy- Casimir  method  (see  Chapter  5). 

7.2.1  Some  Symbolic  Manipulation  Programs 

A  general  purpose  computer  program  for  the  automatic  generation  of  multibody 
systems  in  space,  was  created  using  the  symbolic  manipulation  language  MACSYMA. 
The  program  uses  the  augmented  body  approach  as  given  by  Wittenburg  [59]  to 
symbolically  generate  the  dynamical  equations.  Dynamical  equations  for  a  three 
dimensional  multibody  satellite  connected  in  the  form  of  a  tree  structure  can  be 
generated.  Individual  bodies  are  modeled  as  rigid  bodies,  and  are  connected  to  each 
other  by  revolute  joints  with  as  many  as  six  degrees  of  freedom  per  joint.  Input  to  the 
program  consists  of  the  physical  parameters  of  the  satellite  bodies  like  -  mass,  inertia 
tensor,  relative  locations  of  the  centers  of  mass  and  the  joints;  and  the  position  of 
the  body  1  center  of  mass  in  the  inertial  coordinate  frame.  The  external  torques 
and  forces  acting  on  the  satellite  bodies  could  be  given  either  in  the  inertial  or  the 
body  frame  of  coordinates.  The  output  is  a  set  of  second  order  ordinary  differential 
equations  representing  the  dynamics  of  the  system.  The  dynamical  equations  of  a 
planar  two-body  example  in  Quartararo  [46]  were  generated  using  this  program,  and 
the  results  were  found  to  be  in  complete  agreement. 

Another  symbolic  manipulation  program  to  generate  the  dynamical  equations  of 
multibody  systems  in  space  was  also  implemented.  Here  we  use  a  Newton-Euler 
approach  which  is  based  on  the  method  of  nested  bodies  due  to  Velman  [57].  Details 
of  the  method  are  available  in  Frisch  [13].  This  program  was  also  run  with  the  previous 
example  and  the  results  were  found'to  be  in  total  agreement  with  the  augmented  body 
approach. 
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CHAPTER 

EIGHT 

CONCLUSION 


Dynamics  and  control  of  multibody  systems  has  been  the  focus  of  this  thesis.  It  has 
been  our  contention  that  a  good  understanding  of  the  dynamics  of  multibody  systems 
is  imperative  for  formulating  control  strategies  and  for  the  eventual  goal  of  controlling 
such  complex  systems.  We  have  dealt  with  a  subclass  of  multibody  systems  namely 
planar  multibody  systems  connected  in  the. form  of  an  open  kinematic  chain.  This 
assumption  is  not  restrictive  in  the  sense  that  lot  of  applications  could  be  found. 

A  Hamiltonian  formulation  of  the  dynamics  of  a  general  tree  connected  multibody 
system  was  found.  Existing  translational  and  rotational  symmetries  in  the  problem 
reduce  the  space  on  which  the  dynamics  evolve.  A  Poisson  structure  (bracket)  on 
the  reduced  space  was  found  and  the  dynamics  represented  in  a  elegant  way  using 
these  brackets.  Any  external  and  internal  torque  acting  on  the  system  could  be  taken 
into  account  in  the  modeling.  All  this  was  presented  in  Chapter  three. 

Some  interesting  examples  of  few-body  systems  -  two,  three  and  N-body  (chain), 
were  worked  out  in  Chapter  four.  Equilibria  of  the  multibody  system  dynamics  on  the 
reduced  space  were  studied  in  Chapter  five.  Non-degenerate  equilibria  (rotation  of  all 
the  bodies  in  the  system  with  a  constant  angular  velocities)  were  of  interest.  Stability 
of  such  equilibria  was  examined  using  energy-Casimir  method.  A  formulation  of  the 
energy- Casimir  method  for  a  general  case  of  an  N-body  system  was  given.  It  was 
shown  that  for  the  planar  two-body  example  there  were  two  equilibria  -  a  stretched 
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out  one  (stable  one)  and  a  folded  one  (unstable  one).  A  graphic  computer  simulation 
study  of  the  problem  was  also  presented.  We  have  also  proved  that  there  can  be 
as  many  as  six  equilibria  for  the  special  kinematic  case  of  planar  three-body  system 
where  the  system  center  of  mass  of  the  second  body  is  along  the  line  joining  the  two 
joints.  The  equilibria  of  an  N-body  system  where  all  the  bodies  were  symmetrical 
was  studied  and  it  was  proved  that  the  stretched  out  configuration  is  defintely  an 
equilibrium.  A  numerical  experiment  conducted  to  find  the  stability  of  equilibria 
using  the  energy- Casimir  method  for  specific  cases  of  number  of  bodies  ranging  from 
N  =  2  to  N  =  10  indicated  that  the  stretched  out  configuration  is  stable. 

Some  important  control  problems  were  tackled  in  Chapter  six.  Global  controlla¬ 
bility  of  planar  multibody  systems  was  proved,  given  all  joint  torques  plus  an  external 
torque  (say  a  gas  jet)  on  one  of  the  bodies.  Also  the  multibody  system  could  be  also 
I/O  linearized  by  a  feedback  of  the  form  u(x)  =  a(x)  +  /3(x)v  where  u  is  the  con¬ 
trol  input  of  only  joint  torques.  A  stability  theorem  was  proved  which  states  that  the 
system  converges  to  a  maximal  invariant  subset  consisting  of  all  the  system  equilibria 
on  the  introduction  of  a  feedback  of  the  form  T,  =  C  —  where  C,  >  0. 

Chapter  seven  was  the  concluding  chapter  of  the  thesis  and  presented  our  efforts  in 
applying  symbolic  and  LISP  based  computation  as  tools  for  analysis  and  simulation 
of  planar  multibody  systems.  OOPSS  -  an  Object  Oriented  Planar  System  Sim¬ 
ulator  based  on  object  oriented  programming  methodology,  with  a  general  purpose 
architecture  for  generation,  simulation  and  animation  of  planar  multibody  systems 
was  proposed.  OOPSS  has  been  implemented  for  two  and  three-body  systems  in 
Zeta-lisp  on  a  Symbolics  3600  series  machine  4.  Forcing  inputs  like  the  joint  torques 
could  be  introduced  and  equilibria  could  be  displayed.  A  stabilizing  input  illustrates 
our  stabilization  theorem. 

Future  Work 

We  plan  to  extend  the  above  formalism  to  three-dimensional  systems  i.e.,  for 
systems  with  one  or  more  degrees  of  rotational  freedom  between  two  bodies.  Some 
4 manufactured  by  Symbolics  Inc,  Cambridge,  Mass. 
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work  has  been  done  in  this  regard  by  Grossman,  Krishnaprasad  and  Marsden  [16], 
where  they  consider  two  rigid  bodies  connected  by  a  ball-in  socket  (three  degrees  of 
freedom)  joint.  Equilibria  of  complex  multibody  systems  needs  to  be  identified  and 
stability  issues  studied.  Extension  of  this  modeling  methodology  to  flexible  systems 
are  interesting  to  consider.  One  also  would  like  to  see  a  “controllability”  theorem 
on  the  symplectic  leaf  given  the  joint  torques  only.  Implementation  of  OOPSS  for 
complex  multibody  systems  is  being  considered. 
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Appendix  1 


We  give  the  results  of  a  numerical  experiment  conducted  to  study  the  stability  of 
the  stretched  out  equilibrium  (0»+i,»  =  0,  *  =  1  ,N—  1)  for  the  symmetrical  N-body 
(chain)  example  using  energy- Casimir  method.  The  experiment  was  conducted  with 
the  number  of  bodies  varying  from  2  to  10. 

Finding  the  definiteness  of  the  Hessian  dP[H  +  C)  at  the  above  equilibrium  is 
equivalent  to  finding  whether  the  matrix  given  by  equation  (6.6)  Chapter  5,  is  definite 
or  not  (since  the  psuedo  inertia  matrix  J  is  always  positive  definite).  Call  this  matrix 
M.  Recall  that  the  matrix  M  is  a  (N  —  1  X  N—  l)  matrix.  The  numerical  results  are 
given  below.  Note  that  p  —  d2mu2,  where  d  is  the  distance  between  two  adjacent 
joints,  m  is  the  mass  of  each  body  and  u0  is  the  angular  velocity  with  which  all  the 
bodies  in  the  system  rotate  at  equilibrium. 

N  =  2 

First  principal  minor  of  M  —  0.125. 

So  matrix  M  and  thus  d2(H  +  C )  is  positive  definite. 

N  =  3 

First  principal  minor  of  M  =  0.333p. 

Second  principal  minor  of  M  =  0.10416p2. 

So  matrix  M  and  thus  d2[H  +  C )  is  positive  definite. 

N  =  4 

First  principal  minor  of  M  =  0.5625 p. 

Second  principal  minor  of  M  =  0.5 p2. 

Third  principal  minor  of  M  =  0.25 ps. 
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So  matrix  Af  and  thus  d2{H  +  C )  is  positive  definite. 
N  =  5 

First  principal  minor  of  Af  =  0.8 p. 

Second  principal  minor  of  Af  =  0.12375 p2. 

Third  principal  minor  of  Af  =  0.17875ps. 

Fourth  principal  minor  of  Af  =  0.12289p4. 

So  matrix  Af  and  thus  d2(H  +  C)  is  positive  definite. 
N  =  6 

First  principal  minor  of  M  =  1.041p. 

Second  principal  minor  of  M  =  2.33 p2. 

Third  principal  minor  of  M  =  5.906p8. 

Fourth  principal  minor  of  Af  =  11.813p4. 

Fifth  principal  minor  of  Af  =  10.336/?5. 

So  matrix  M  and  thus  d2(H  +  C)  is  positive  definite. 
N  =  7 

First  principal  minor  of  M  =  1.286 p. 

Second  principal  minor  of  Af  =  3.795p2. 

Third  principal  minor  of  Af  =  13.964ps. 

Fourth  principal  minor  of  Af  =  49.093/?4. 

Fifth  principal  minor  of  Af  =  125. 46p5. 
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Sixth  principal  minor  of  M  =  133. 3p6. 

So  matrix  M  and  thus  d2(H  +  C )  is  positive  definite. 

N  =  8 

First  principal  minor  of  M  =  1.531  p. 

Second  principal  minor  of  M  =  5.625p2. 

Third  principal  minor  of  M  —  27.34/j3. 

Fourth  principal  minor  of  M  =  140.0p4. 

Fifth  principal  minor  of  M  =  630.0ps. 

Sixth  principal  minor  of  M  =  1960.0p6. 

Seventh  principal  minor  of  M  =  2450.0p7. 

So  matrix  M  and  thus  d2(H  +  C )  is  positive  definite. 
N  =  9 

First  principal  minor  of  M  =  1.777 p. 

Second  principal  minor  of  M  =  7.826 p2. 

Third  principal  minor  of  M  =  47.44ps. 

Fourth  principal  minor  of  M  =  321. 19p4. 

Fifth  principal  minor  of  M  =  2107 p5. 

Sixth  principal  minor  of  M  =  11555p6. 

Seventh  principal  minor  of  M  =  42370 p7. 

Eighth  principal  minor  of  AT  =  60908p8. 
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So  matrix  M  and  thus  d?(H  -f  C)  is  positive  definite. 


N  =  10 

First  principal  minor  of  M 

=  2.025 p. 

Second  principal  minor  of  M 

=  10.4p*. 

Third  principal  minor  of  M 

=  75.64ps. 

Fourth  principal  minor  of  M 

=  639. lip4. 

Fifth  principal  minor  of  M 

=  5548 p5. 

Sixth  principal  minor  of  M 

=  44382p6. 

Seventh  principal  minor  of  M 

=  287101p7. 

Eighth  principal  minor  of  M 

=  1212205p8. 

Ninth  principal  minor  of  M 

=  1969834p®. 

So  matrix  M  and  thus  d2(H  +  C )  is  positive  definite. 
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Appendix  2 


Col  1 


Col  2 


MMNMHMMHNMNMMMMMHMHMNMMHHKMMMMttttHMMMMHMMNI 

DYNAMICS  OF  PLANAR  MULTIBODY  SYSTEMS 

NNKNNMNIfNKNHIINMNtlMMItlllfNNintHMMKNNNMIlMMMMlIHI 


Ths  Multibody  system  data  la  In  flla  < 
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The  path  fron  center  of  nasa  of  body  1  to  body  2  la  Cl,  2] 
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The  J  matrix  Is  given  below 
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H  M  C06( THETA  .  THETB  >  C  1  B  1 
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